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19.  Abstract 

By  mapping  the  ray  paths  and  the  reflection  points  of  rays  which  set  out 
with  specified  phase  velocities  a  semiquantitative  measure  of  the 
disruption  of  Lg  propagation  patterns  by  complex  structure  can  be  made.  In 
addition  it  is  possible  to  estimate  the  likely  extent  of  conversion  between 
Rayleigh  and  Love  wave  propagation.  This  method  has  been  successfully 
applied  to  predicting  patterns  of  Lg  propagation  in  central  Asia  and  the 
southwestern  United  States. 

A  recent  sequence  of  earthquakes  in  the  Northern  Territory  of  Australia 
provided  an  opportunity  to  record  regional  seismic  phases  in  a  Pre  Cambrian 
geologic  environment  (similar  to  a  number  of  areas  within  the  USSR) . 

The  seismograms  display  a  high  Lg  velocity  (around  3.7  km/s)  and  the 
attenuation  of  the  regional  phases  with  distance  is  incur. sistent  with  the 
geometric  spreading  relations  commonly  adopted  for  paths  in  the  United 
States  and  Europe.  This  discrepancy  appears  to  arise  fm.m  the  gradational 
character  of  the  crust  -  mantle  transition  in  central  Australia. 

It  is  therefore  importa.it  that  corrections  tor  attenuation  taxe  account 
of  the  character  of  the  crustal  structure  because  inappropriate  corrections 
could  build  in  systematic  error  to  estimates  of  source  magnitude  (and  hence 
yields)  derived  from  regional  observations. 


Accession  For 

"ntis  gra&i  c 

DTI(3  TAB  ^ 

Una-.uiounced  | 

.lut  tif  Icntlon _ 

□ 

E"’ _ 

?■  Ir;:  ,:i:i  '.n/ _ _ 

"1.1  •  ■ "  i  ./  Codes 


±L 


Rl^GlONAL  PHASE  PROPAGATION  IN  THREE  DIMENSIONAlLY  HETEROGENEOUS 
Media 


RESEARCH  OBJECTIVES: 

The  aim  of  this  work  is  to  develop  techniques  which  can  be  used  to  describe  the 
propagation  of  regional  seismic  phases  in  three  dimensionally  varying  media,  in 
order  to  improve  understanding  of  the  nature  of  such  phases  and  the  way  in 
which  the  characteristics  of  the  seismic  source  can  be  modified  by  propagation  to 
the  receiver.  Such  information  will  be  valuable  in  assessing  the  effects  of 
geological  structure  on  the  behaviour  of  potential  discriminants  between  nuclear 
explosions  and  earthquakes. 

The  problem  of  discrimination  becomes  more  difficult  for  small  events  recorded 
at  considerable  distances.  The  regional  seismic  phases,  from  which  the 
discriminant  properties  are  extracted,  will  have  passed  through  regions  of 
considerable  horizontal  variability  in  seismic  properties  in  their  passage  from 
source  to  receiver. 

VVe  have  therefore  set  out  to  find  ways  to  improve  the  description  of  the 
interaction  of  regional  phases  with  crustal  structure  via  a  combination  of 
theoretical,  computational  and  observational  work. 

research  status  and  ACHIEVEMENTS: 

The  crustal  waveguide  is  a  zone  of  complex  variations  in  seismic  wavespeeds. 
Although  it  is  feasible  to  generate  specific  models  of  the  heterogeneity  structure 
it  is  more  economical  to  use  a  stochastic  description  of  the  structure.  However, 
with  a  stochastic  model  it  is  not  possible  to  use  the  coupled  modes  techniques 
developed  by  Kennett  (1984)  for  the  description  of  regional  S  waves  (the  phases 
Lg,  Sn).  Fortunately  it  is  possible  to  modify  techniques  originally  devised  for 
wave  propagation  in  optical  fibres  to  generate  a  new  method  for  estimating  the 
structural  component  of  the  apparent  attenuation  of  guided  waves  with  distance. 

For  a  particular  stochastic  model  of  crustal  and  mantle  heterogeneity  it  is  now 
possible  to  calculate  the  attenuation  of  different  modes  as  a  function  of  distance, 
including  intrinsic  attenuation  and  scattering  loss  by  transfer  to  other  modes  of 
propagation.  Such  statistical  estimates  of  attenuation  can  then  be  compared  with 
observations  of  the  attenuation  of  regional  phases  for  many  paths  traversing  a 
single  region.  For  Lg  and  Sn  scattering  attenuation  can  represent  nearly  a  third 
of  the  apparent  attenuation. 

A  full  description  of  the  coupled  mode  technique  for  stochastic  media  and  the 
resulting  estimates  of  scattering  attenuation  are  given  in  the  attached  paper: 
Guided  wave  attenuation  in  laterally  varying  media. 
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The  Lg  phase  is  one  of  the  most  prominent  features  of  many  regional 
seismograms  in  the  frequency  band  from  0.5  to  3.5  Hz  and  its  size  makes  it 
valuable  for  the  detection  of  small  events.  However  the  amplitude  of  Lg  is 
known  to  be  sensitive  to  changes  in  crustal  structure  and  this  behaviour  reduces 
its  utility  for  discrimination  purposes.  Given  a  particular  source  location  it  is 
therefore  highly  desirable  to  be  able  to  predict  the  effects  of  crustal  heterogeneity 
on  the  characteristics  of  Lg  propagation.  A  large  part  of  the  effect  of  3-D  structure 
arises  in  the  shape  of  the  crustal  waveguide  and  by  using  a  ray  theoretical 
description  of  the  Lg  phases,  topography  at  the  surface  and  the  Moho  can  be 
included.  Rays  are  traced  at  different  initial  angles  to  the  vertical  depending  on 
the  phase  velocity  under  consideration.  The  rays  are  monitored  in  plan  view 
and  by  mapping  the  reflection  points  at  the  surface  and  Moho,  a 
semiquantitative  measure  of  the  disruption  of  Lg  can  be  made.  Also,  the  tilt  of 
the  rays  can  be  used  to  judge  the  degree  of  interconversion  between  Love  and 
Rayleigh  wave  energy. 

The  development  of  the  technique,  testing  on  simple  models  and  application  to 
the  complex  structures  in  central  Asia  are  described  in  the  accompanying  paper; 
The  effect  of  three  dimensional  structure  on  Ly  propagation  patterns. 

Despite  the  relative  simplicity  of  the  method  the  major  characteristics  of  Lg 
propagation  are  well  predicted.  By  including  the  possibility  of  secondary  sources 
scattering  due  to  specific  structural  features  can  be  included  and  this  approach  is 
exploited  in  the  section; 

Guided  wave  tracking  in  3-D  -  a  tool  for  interpreting  complex  regional 
seismograms 

in  which  the  Lg  wave  tracking  procedure  is  applied  to  the  southwestern  United 
States.  Predicted  strong  scattering  from  the  Sierra  Nevada,  reflection  from  the 
ocean-continent  transition  and  resonance  within  the  narrow  corridor  of  Baja 
California  match  up  well  with  features  on  observed  seismograms. 

Available  detailed  observations  of  regional  seismic  phases  cover  a  relatively 
limited  range  of  tectonic  environments  and  there  is  therefore  the  possibility  that 
many  well  accepted  observations  may  be  specific  to  certain  styles  of  crustal 
structure.  The  Tennant  Creek  earthquake  sequence  in  the  Northern  Territory  of 
Australia  in  an  area  of  PreCambrian  outcrop  provided  an  opportunity  to 
investigate  regional  phase  propagation  in  a  different  context.  Recordings  were 
made  on  a  roughly  linear  area  along  the  highway  running  between  the  arrays  at 
Warramunga  and  Alice  Springs  giving  sparse  but  useful  coverage  over  a  400  km 
distance  span. 

The  apparent  velocity  of  Lg  is  quite  fast  (3.7  km/s)  and  its  decay  pattern  is  not 
consistent  with  the  conventional  assumption  of  geometric  spreading  in  the  form 
r-5/6  These  unusual  features  appear  to  be  related  to  the  gradational  character  of 
the  crust-mantle  transition  in  this  region.  The  accompanying  report 
The  effects  of  crustal  velocity  gradients  on  the  propagation  of  regional  S  phases 
describes  the  analysis  of  this  data.  Work  continues  on  attempting  to  characterise 
the  dependance  of  the  parameters  of  the  regional  phases  in  the  nature  of  the 
crust-mantle  tiansition. 
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Guided  Wave  Attenuation  in  Laterally  Varying  Media 


B. L.N'.  Kennett 

RcMirch  School  of  Earth  Sciences,  Australian  National  University, 

C. P.O.  Box  4,  Canberra  ACT  2601,  Australui 


Summary 

Coupled  mode  techniques  for  guided  wave  propagation  are  extended  to  two- 
dimensional  stochastic  heterogeneity  superimposed  on  a  stratified  medium.  This 
approach  requires  the  variations  to  be  smoothly  varying  and  of  modest  size  (less  than 
±2  per  cent).  By  averaging  over  an  ensemble  of  statistically  similar  models  coupled 
equations  for  the  modal  energy  transport  can  be  generated.  The  intermode  coupling 
depends  on  the  horizontal  correlation  functions  for  the  heterogeneity  in  the  crust  and 
mantle,  and  the  integrated  effect  of  the  vertical  variations  in  velocity  and  the  modal 
eigenfunctions. 

For  a  particular  stochastic  model  the  attenuation  of  a  single  mode  as  a  function  of 
distance  can  be  calculated  as  a  superposition  of  intrinsic  attenuation  and  scattering 
loss  by  energy  transfer  to  other  modes  of  propagation.  These  statistical  estimates  of 
attenuation  can  be  compared  with  observations  of  regional  phases  travelling  over  a 
variety  of  paths  in  a  single  region.  ITir  Lg  and  Sn  phases,  intermode  scattering  may 
represent  up  to  30  per  cent  of  the  apparent  loss. 

Ka/  Words:  Attenuation,  Guided  wave  propagation,  stochastic  heterogeneity,  regional 
seismic  phases 


Introduction 

Detailed  studies  of  the  structure  of  the  crust  and  upper  mantle  show  lateral 
heterogeneity  on  a  wide  variety  of  scales  superimposed  on  a  basic  stratification  with 
depth.  Guided  seismic  waves  travelling  nearly  horizontally  through  such  structures 
are  particularly  vunerable  to  scattering  due  to  the  presence  of  heterogeneity. 

A  convenient  way  of  representing  such  guided  waves  is  in  terms  of  a 
superposition  of  surface  wave  modes  whose  eigenfunctions  are  largely  confined  to  the 
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waveguide.  Such  a  description  has  been  used  by  Malin  (1980),  Wang  &  Herrmann 
(1988)  in  simulations  of  the  coda  of  local  earthquakes.  For  longer  distance 
propagation  Malichewsky  (1987)  has  summarised  the  techniques  for  handling  the 
interactions  of  surface  waves  with  vertical  interfaces  or  sharp  transitions.  These 
methods  relate  the  amplitudes  of  the  surface  wave  modes  in  the  structure  on  one  side 
of  the  boundary  to  those  on  the  other  side.  In  a  similar  way,  for  continuous 
heterogeneity  the  wavefield  can  he  described  in  terms  of  surface  wave  modes  if 
allowance  is  made  for  energy  transfer  between  modes  (Kennett  1984;  Maupin  1988). 
These  coupled  mode  methods  have  generally  been  fonnulated  for  deterministic 
velocity  structures  and  much  less  work  has  been  done  on  a  stochastic  description  of 
guided  wave  propagation  even  though  this  has  been  extensively  developed  for  body 
waves  (see  e.g.  Aki  &t  Richards  1980,  chapter  13;  Hudson  1982). 

We  here  present  an  adaptation  of  the  coupled  mode  procedure  to  a  2-D 
stochastic  heterogeneous  medium,  and  generate  coupled  t‘<.]uations  for  the  modal 
energy  transport  averaged  over  an  ensemble  of  statistically  similar  models.  The 
approach  is  based  on  techniques  developed  for  fibre  optics  (Marcuse  1974),  but  is 
extended  to  multiple  waveguides  e.g.  the  crust  and  upper  mantle  for  regional  seismic 
phases.  For  a  specified  stochastic  model  of  the  velocity  heterogeneity  we  are  able  to 
estimate  for  each  mode  the  contribution  to  the  modal  loss  factor  iQ„((o)  induced  by 
scattering  into  other  modes  of  propagation.  Such  estimates  can  be  compared  with 
observations  of  guided  waves  over  a  variety  of  paths  in  a  single  geographic  region,  as 
e.g.  the  observations  of  Nuttli  (1980)  for  Sn  and  Lg  phases  propagating  across  Iran. 

Coupled  Mode  Equations: 

Kennett  (1984)  has  shown  how  the  dispbcement  and  traction  fields  for  guided 
seismic  waves  in  2-D  laterally  varying  media  can  be  represented  as  a  sum  of  modal 
eigenfunctions  with  horizontally  varying  coefticients.  When  ihe  variation.^,  in  seismic 
parameters  do  not  show  a  systematic  trend  horizontally,  the  expansion  may  be  made 
in  terms  of  the  modes  for  a  fixed  reference  structure.  Maupin  &  Kennett  (1987)  have 
shown  how  the  effect  of  inclined  interfaces  can  be  included  in  the  .scheme  bv 
modifying  the  matrix  elements  in  the  differential  equations  for  coupling  between  the 
modal  coefficients.  When  the  structure  has  systematic  horizontal  variations  a  more 
effective  representation  is  to  work  in  tenns  of  the  local  modes  of  the  staicture 
(Maupin  1987).  Although  the  nature  of  the  coupling  coefficients  differ,  the  two  styles 
of  modal  representation  lead  to  sets  of  coupled  equations  with  similar  structure. 
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The  coupled  mode  equations  give  a  representation  in  terms  of  a  particular 
velocity  structure.  However,  in  considering  different  classes  of  heterogeneity  it  can  be 
advantageous  to  adopt  a  stochastic  treatment  for  the  heterogeneity  structure.  We  will 
consider  the  case  where  the  heterogeneity  has  no  more  than  ±2  per  cent  variation  in 
seismic  parameters  from  those  of  a  reference  model.  This  provides  a  good 
repre.sentation  of  many  situations  in  the  crust  and  uppermost  mantle  and  can  be  well 
treated  using  the  coupled  mode  approach  with  a  fixed  set  of  reference  modes  (Kennett 
19S9). 

At  each  point  in  the  2-D  varying  medium  we  represent  the  displacement  field 
w(x,z)  as  a  sum  of  vertically  varying  modal  eigenfunctions  with  coefficients  which 
depend  only  on  horizontal  position  (Kennett  1984).  Thus 

w(,r ,:  )  =  5^  c,  (.r  )  exp  { i^  r  )  ,2),  ( 1 ) 

r 

where  w“  is  the  eigenfunction  for  the  rth  mode  with  horizontal  wavenumber  k^,  and 
the  sum  is  taken  over  both  forward  and  backward  travelling  modes.  The  total 
number  of  modes  must  be  chosen  large  enough  to  include  all  the  propagation 
processes  of  interest  (see  Maupin  &  Kennett  1987).  The  eigenfunctions  for  the  modes 
are  normalised  so  that  they  have  the  same  horizontal  energy  transport. 

The  horizontal  evolution  of  the  modal  coefficients  is  described  by  a  set  of 
ordinary  differential  equations 

-:^Cr(.r)  =  X  Kri{x)e\p[i%-kr)x  )Cs(a-),  (2) 

ox 

where  the  coupling  matrix  K  depends  on  the  deviations  of  the  properties  of  the  actual 
medium  from  those  of  the  reference.  Using  the  concise  notation  of  Maupin  &  Kennett 
(1987),  which  is  suitable  for  anisotropic  media 

Kf5  =  i  j ciz  +  Apw'w^.w^ 


-  w^.AQiidj  wf  +  w°.A(C3]C]]  )1}'s  -  t“r.A(Cii  C  n)d^  Ws°] 


+  (3) 

n 

where  tf,  is  the  horizontal  traction  derived  from  w°  and  we  have  w’ritten 
w"  =  w°(-4  ,z ).  In  terms  of  the  elastic  modulus  tensor 

(C,f)ki  =Ck,)j,  Q:\:<  = 
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The  intci  *  .ce  temis  depend  on  the  slope  of  the  interface  and  the  jumps  in  the 
horizontal  traction  1“^  across  the  interface. 

The  set  of  equations  (2)  include  the  possibility  of  both  reflection  and 
transmission.  However,  dirc'Ct  calcubtions  using  the  coupled  mode  ec^uations  show 
that  reflected  waves  can  be  neglecttxl  without  appreciable  error,  provided  that 
deviations  in  the  seismic  parameters  from  the  reference  model  are  not  too  large  (less 
than  ±2  per  cent)  and  the  scale  of  variation  is  not  rapid  compared  with  the  horizontal 
wa\elengths  involved.  This  neglect  of  reflected  waves  can  also  be  justifietl  bv  using 
first-order  Born  scattering  results  in  the  wavenumber  domain  (Kennelt  1972).  I'or  a 
heterogeneity  pattern  with  wavenumber  spectrum  /  (A),  the  scattering  between  two 
wavenumbers  and  is  proportional  to  I'or  reflected  waves  k^,  is 

negative  and  so  the  differetice  k^-kf,  will  be  large.  Thus  scattering  into  reflected 
waves  will  onlv  be  noticeable  if  f  ik)  has  significant  amplitude  for  large  k,  i.e.  if  the 
heterogeneity  itself  varie^  on  very  small  horizontal  scales  (or  is  di.scontinuous). 

We  will  therefore  restrict  attention  to  transmitted  waves  and  so  consider  the 
limited  set  of  coupled  equations 

-r-Cr{x)='£  K,H.t)exp|i(A;-Av).v  )Cj{.r).  (4) 

With  the  restrictions  we  have  placed  on  the  nature  of  the  heterogeneity  a  good 
approximation  for  the  modal  coefticients  after  passage  through  a  small  horizontal 
distance  h  is  given  by 

Lrix+h)  =  Crix)  +  Y.  ‘slv)  j  dq  K,, (rj)exp(i(A, -A,)(j  ),  (3) 

-  ti  t 

the  tirst  temi  in  a  systematic  expansion  in  terms  of  the  coupling  coefficients 
higher  order  interactions  can  be  neglected  since  we  have  assumed  weak  heterogeneitv. 
We  are  thus  able  to  treat  the  set  (cy)  as  the  amplitude  distribution  across  the  modes 
at  the  horizontal  position  x. 

Coupled  Power  Equations 

The  set  of  coupled  mode  equations  (2)  can  be  solved  for  a  particular 
heterogeneity  model  and  then  provides  a  detailed  description  of  the  phase  and 
amplitude  behaviour  of  all  the  modes  at  each  point  .r  in  the  waveguide.  Bv 
combining  this  information  for  many  fre(|uencies  we  can,  in  principal,  constnict 
theoretical  seismograms  but  very  substantia!  computational  effort  is  ri\]uired 
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We  can,  however,  get  a  measure  of  the  energy  redistribution  between  modes  as  a 
function  of  distance  in  propagating  through  a  random  heterogeneity  by  considering 
the  average  power  over  an  ensemble  of  heterogeneity  models.  The  resulting 
equations  for  the  modal  power  contributions  for  each  frcKjuency  are  a  coupled  set  of 
first  order  differential  equations  with  constant  and  symmetric  coefficients. 

We  will  adopt  a  physically  based  approach  to  the  derivation  of  the  coupled 
power  equations  following  Marcuse  (1974).  The  results  we  derive  can  be  justified 
formally  using  the  techniques  of  stochastic  differential  equations  described  in  Kohler 
&  Papanicolaou  (1977). 

We  have  adopted  a  normalisation  for  modal  eigenfunctions  in  which  the 
horizontal  energy  transport  is  equal  for  each  mode.  As  a  result,  for  a  particular 
heterogeneity  model  a  measure  of  the  power  in  an  individual  mode  P„  is 


and  from  (4)  satisfies  the  horizontal  evolution  equation 


()  n  '  ^ 

—  P  =  r  —  c 

dx  '"ax  " 


X  •  ^  ■ 


=  C^'L  Kmn exp { i(l:„ -k„ )x  }  +  c.c. ,  (7) 

n 

wliere  c.c.  denotes  the  complex  conjugate  of  the  previous  term.  We  now  consider  an 
ensemble  average  over  a  collection  of  statistically  similar  heterogeneity  models,  built 
according  to  the  same  prescription  but  not  identical.  In  particular  we  assume  that  the 
phases  of  any  periodic  components  in  the  heterogeneity  are  randomly  distributed 
across  the  members  of  the  statistical  ensemble. 

The  averaged  power  in  each  mode 


and  its  change  with  horizontal  position  will  be  governed  by  the  ensemble  average  of 
equation  (7): 


_a 

dx 


T<K  mn  c„>cxp{i(*:„-A„,)x)  +  c.c. 


(9) 


We  will  assume  that  each  of  the  elements  defining  the  heterogeneity  Ap,  AC  n  , 
■\Qw,  A(Cf,'Cn),  A(C3]Cn  )  have  similar  statistical  properties  and  the  heterogeneity 
values  at  widely  separated  points  are  uncorrelated.  The  coupling  coefficients  will  be 
described  by  a  stationary  random  process  with  a  finite  correlation  length  D  in  the 
horizontal  direction.  We  will  assume  the  average  heterogeneity  level  across  the 
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ensemble  vanishes  so  that 


<K„„  (x  )>  =  0. 

For  the  guided  waves  travelling  in  the  direction  of  increasing  x  wc  anticipate 
that  the  modal  amplitude  c„(x)  and  the  coupling  matrix  K  (x)  will  be  uncorrelated  if 

x'-x  3>  D, 

and  then 

<c„  (x'  k„  (x'  )K„„  (x  )>  =  <c„  (x'  )c„  (x  (x  )> 

with  similar  behaviour  for  products  of  coupling  coefficients.  We  can  exploit  this  result 
in  (9)  if  we  express  the  modal  field  at  x  in  terms  of  that  at  x'  using  (5).  Under  the 
assumption  that  third  order  interaction  terms  of  the  type  <KKK>  can  be  neglected 
because  the  heterogeneity  is  small,  we  find 

=  £[  <c^(x')c,{x')>e\p{i(k„-k„)x)  jdq  <K„„(x)K„,{q)>e\p{\(k,-k„)x  ] 

n,r  X 

x' 

+  <c„(x'  kr(x'  )>exp{i(k„-k„)xjjdq  <K„„{x)K^r^q)> exp[-i(kr-k„)x  ) 

X 

+  C.C.].  (10) 

Since  we  have  assumed  a  short  correlation  length  we  can  extend  the  lower  limit  of 
integration  to  and  so  recast  the  integrals  in  the  form  of  half  range  Fourier 
transforms  over  the  correlation  between  the  coupling  coefficients 

x' 

exp[i(k„-k„)x]  jdq  <K„„{x)K^,{q)> exp ( i(*:„, -k, )q  ) 

X 

oo 

=  expli(;c„-<:,).r  }|d»  <K„„(x)K^,(x-M)>cxp(-i(fc„-A.y)u  }.  (11) 

The  coupled  power  equations  can  therefore  be  written  as 

<cJl(.r')t%(x')>exp{i(fc^-<:„)x)|dH  <K„„(x)K„,{x-u)>exp[\(k„-k,)u 

+  <c„  (x'  )c/(.v')>exp(i(lc„-A,)x  )  |du  <K„„(x  )K^,{x-u  )>exp(i(i.,-J:„  )h 

+  C.C.].  (12) 
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We  anticipate  that  the  modal  amplitudes  will  tend  to  have  random  phase  so  that 
<c„c^>  =  <\c„  12>5„„ 

and  this  selection  of  diagonal  elements  will  be  reinforced  by  a  situation  akin  to 
stationary  phase:  the  main  contribution  from  the  right  hand  side  of  (12)  will  arise 
from  non-oscillatory  terms  in  the  integrals.  Hence  we  are  able  to  drop  one 
summation  to  give 

oo 

=  zt  s„(x')  fdu  <K„„{x)K„„(x-u)> e\p[i{k„-k„)u  } 

n  n 

oo 

+  S„  (x')  Jdw  <K„„(x)K^„{x-u)>exp{\{k„-k„)u  ) 


We  have  assumed  that  the  modal  coefficients  are  slowly  varying  in  going  from  .r  to 
x' ,  and  so  S„  will  be  nearly  constant  over  this  interval.  As  a  result,  we  can  replace 
S„(x')  by  Sm(x)  on  the  right  hand  side  of  (13).  Further,  we  require  that  the  total 
power  should  be  independent  of  x  in  a  lossless  medium  so 

k„„(x)  =  -k;„(x). 

Thus 

^S„,  =  jdu  {<K„,„(x)KJ;„(x-h)>(S„(x)-S„(x)])  exp[i(l:„-lr,„)«], 

where  we  have  used  the  fact  that  the  term  in  braces  is  real  to  incorporate  the  complex 
conjugate  term  as  an  integral  along  the  negative  «  axis.  (Note  that  in  this  perfectly 
elastic  situation  after  propagation  through  a  very  large  distance  the  power  in  all  the 
modes  will  equalise  to  the  same  level.) 

Equation  (14)  covers  the  general  case  of  weak  heterogeneity  but  we  can  simplify 
further  if  we  make  specific  assumptions  about  the  form  of  the  heterogeneity  model. 
The  near  stratification  within  the  Earth  will  impose  different  stochastic  properties  on 
the  horizontal  and  vertical  variations  of  the  heterogeneity.  The  coupling  coefficients 
K„,„  are  defined  in  terms  of  an  integral  over  the  full  depth  of  the  heterogeneity.  After 
taking  the  ensemble  average  over  the  set  of  statistically  similar  models,  the  evolution 
equations  for  the  average  modal  power  will  be  dominated  by  the  nature  of  the 
variations  in  the  heterogeneity  in  the  horizontal  direction. 
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Separation  of  crust  and  mantle  heterogeneity  contributions 

In  general  we  expect  the  heterogeneity  in  the  crust  to  have  a  different  character 
to  that  in  the  mantle.  Since  we  have  assumed  weak  heterogeneity,  we  are  able  to 
extract  a  common  functional  dependence  from  each  of  the  coupling  coefficients  K„„ 
for  the  horizontal  behaviour  of  the  heterogeneity  in  each  of  the  crust  and  mantle 
zones 

K„„(x)  =  f^  (x  )K^„  (X )  + f'^{x)K!^„(x).  (15) 

is  the  horizontal  variation  function  for  the  crustal  heterogeneity  and  K^„(x)  is 
an  integral  restricted  to  the  crust.  Similarly  f'^(x)  is  the  horizontal  variation  function 
for  the  mantle  and  K/!f„(x)  is  defined  by  an  integral  over  the  span  of  the  mantle.  We 
do  not  require  a  complete  separation  of  the  horizontal  and  vertical  dependence  of  the 
velocity  and  density  perturbations  from  the  reference.  The  part  of  the  variation 
unaccounted  for  by  f‘^(x)  remains  in  the  crust  and  mantle  coupling 

coefficients  KI^„(x)  which  depend  on  the  properties  as  a  function  of  depth  at 

the  location  .r .  When  we  take  the  average  over  an  ensemble  of  statistically  equivalent 
models,  we  remove  this  dependence  on  location.  The  ensemble  average 

<K„„(x)K;„(x-u)>  =  <f^(x)f^(x-u)>H^^J^  +  <f^(x)f^(x-u)>H^j;fj^ 

=  R^(«)lie^j2  +  RW(M)lCn|2,  (16) 

in  terms  of  the  autocorrelation  functions  of  the  horizontal  variation  of  the 

crust  and  mantle  heterogeneity  e.g.  R^ {u)  =  <f^ {x)f^- {x-u)>.  represent 

ensemble  averages  over  the  vertical  variations.  We  .have  made  the  not  unreasonable 
assumption  that  the  horizontal  variations  in  the  crust  and  mantle  heterogeneity  are 
uncorrelaled. 

With  the  expression  (16)  for  the  heterogeneity  contribution  the  evolution 
equations  for  the  averaged  modal  power  can  be  recast  as 

+  iK;!fj^<lF^(k„-k„)l‘>j  (17) 

where  <IF‘'  P>,  <IF^I‘>  are  the  power  spectra  of  the  cnist  and  mantle  correlation 
functions  e.g. 

OB 

<IF^(k„-k„)l^>  =  J  du  R'- (i/)exp{j(k„-k„)u  ).  (18) 
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We  have  so  far  ignored  the  possibility  of  intrinsic  attenuation  of  the  modes  as 
they  propagate  through  the  medium,  but  it  is  easy  to  compensate  for  the  effect  of  Q 
by  introducing  a  simple  power  loss  term  into  (17).  Thus  for  an  attenuative  medium 

=  -k„Q-^S„{x)^  X  H„„[S„ix)-S„{x)],  (19) 

dx  „=o 

where  the  power  coupling  coefficients 

=  \K^„\^<\F<=(k„-k„)\^>  +  \K!:!„\^<\F^ik„-k„)\^>  (20) 

and  is  the  intrinsic  loss  factor  for  the  mth  mode.  The  spatial  loss  factor  is 
related  to  the  temporal  loss  factor  '  for  the  same  mode  by  the  ratio  of  the  group 
velocity  (U„)  and  phase  velocity  (c„)  for  the  mode 

'Q„^'  =  (U„/c„)Q„'. 

The  set  of  coupled  equations  (19)  can  be  readily  solved  numerically  for  any  given 
power  distribution.  Similar  equations  can  be  derived  for  other  configurations  of 
heterogeneity,  the  differences  will  lie  in  the  form  of  the  power  coupling  coefficients  as 
e  g.  where  the  mantle  has  to  be  subdivided  into  different  zones. 

For  small  perturbations  of  an  interface 

K„„{x)  =  f‘{x)KL(x),  (21) 

where  f'{x)  is  the  shape  of  the  varying  interface.  The  coupling  integral 
will  be  confined  to  a  depth  band  representing  the  extent  of  the  interface  variations, 
and  will  depend  on  the  contrast  in  seismic  properties  across  the  interface.  In  this  case, 
the  equivalent  to  (16)  is 

<K„„(.x)K^„{x-u)>  =  <f‘{x)f'{x-u)>\Rj„„\^  =  K'(u  )IK^„I^  (22) 

where  R'  is  the  autocorrelation  function  of  the  interface  shape  and  the  ensemble 
average  over  the  coupling  terms.  The  subsequent  analysis  will  parallel  equations  (17- 
19)  and  the  corresponding  power  coupling  coefficients 

=  \l^lj^<\F‘{k„-k„)C->.  (23) 

For  interfaces  with  considerable  contrast,  e  g.  the  crust-  mantle  boundary  the  effect  of 
random  topography  can  be  comparable  to  that  of  velocity  perturbations.  The  present 
approach  is  not  suitable  for  situations  with  a  systematic  variation  in  the  depth  of  an 
interface. 
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Attenuation  via  intennode  coupling 

We  can  rewrite  the  coupled  power  equations  (19)  in  a  form  which  emphasises 
the  behaviour  of  the  mth  mode 

)S„  (X  )  +  S 'Hmn S„  (X ),  (24) 

n=0 

where  I'  denotes  a  summation  omitting  m=n,  and 

K  =  I'  H„„ 

n  =0 

i.e.  a  sum  over  one  column  of  the  matrix  of  power  coupling  coupling  coefficients  H 
excluding  the  diagonal  element.  Note  that  coupling  a  mode  to  itself  will  not  affect  the 
power  distribution  with  the  natural  consequence  of  excluding  the  diagonal  eLmcn^. 

For  a  given  initial  power  distribution  an  equivalent  integral  equation  to 

(24)  is 

S„(x)  =  e\p[-(k„Q-'+b„)(x-XQ)]So„{xo) 

N 

+  S '  jdi-  exp{-{k„  Q-'+b„  )(x -i; ) )  H„„  S„  (v ).  (25) 

n=0  X(, 


Since  we  have  assumed  relatively  weak  heterogeneity  the  dominant  behaviour  can  be 
seen  from  the  first-order  approximate  solution 


S„(x}=  exp  { -a„  (x  -x  o) )  S (x (,) 
N 

I 

n=fl  X(, 


-t-  X '  j d I'  exp  ( (x  -y ) }  H„„  exp ( -(j„  ( i-  -x  „) )  S„  (x o) . 


(26) 


where  we  have  introduced  the  effective  decay  rate  with  distance  for  energy  in  the  mth 
mode 

=  k„  Q-^  +b„=k„  Q„ '  +  .  (27) 


If  we  now  consider  an  initial  distribution  Son(xo)  concentrated  solely  in  the  mth 
mode,  we  see  from  (25,26)  that  the  main  part  of  the  behaviour  of  the  mth  mode  with 
distance  will  be  decay  as  expl-(j„  (x-Xo)).  There  will  be  limited  retransfer  of  energy 
to  the  mth  mode  by  secondary  scattering  from  other  modes. 

We  can  therefore  characterise  the  energy  decay  rate  for  the  mth  mode  by  a„  and 
regard  the  intermode  coupling  as  introducing  a  'scattering'  attenuation  term  for  the 
Hith  mode 
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(28) 


iQm’  =  k„'b„ 


where  we  have  reinstated  the  explicit  form  (20)  for  the  power  coupling  coefficients 
between  the  modes.  Equation  (26)  represents  the  first  term  in  an  iterative  solution  in 
terms  of  H  and  since  all  elements  of  H„„  are  positive  the  full  solution  will  give  a 
slightly  larger  value  for  S„  (x )  than  estimated  from  (26).  The  loss  factor  j  Q~^  defined 
by  (28)  is  therefore  a  little  too  small,  but  gives  a  very  useful  indication  of  the 
behaviour.  Equation  (28)  represents  the  spatial  loss  factor;  the  corresponding  temporal 
loss  factor  for  the  wth  mode  will  be 

For  regional  seismic  phases,  the  scattering  attenuation  will  be  significantly 
different  for  the  mantle  phases  than  the  crustally  guided  waves.  The  dominant  term 
in  (25)  for  Lg  type  modes  will  be  )l^>  and  since  we  anticipate  the 

largest  heterogeneity  to  be  concentrated  in  the  crust  j  Qm'  can  be  significant.  For  the 
Sn  modes  the  main  contribution  wilt  be  from  the  mantle  heterogeneity. 

When  we  wish  to  compare  the  behaviour  of  a  set  of  modes  at  fixed  frequency  o), 
it  can  be  advantageous  to  work  with  the  energy  decay  factors  a„,  b„  since  these  terms 
are  directly  comparable  between  modes. 

The  expressions  for  the  scattering  attenuation  for  individual  modes  derived  in 
this  section  are  dependent  on  our  assumption  of  2-D  heterogeneity.  For  3-D 
heterogeneity  there  is  the  additional  complication  that  there  is  the  likelihood  of  energy 
being  scattered  in  directions  other  than  forwards  and  backwards  along  the  local 
direction  of  propagation.  Once  scattered  out  of  a  particular  mode,  energy  is  unlikely 
to  return  to  it  in  the  case  of  3-D  heterogeneity;  thus  attenuation  will  be  stronger  than 
the  estimates  based  on  our  2-D  model. 

However,  when  we  consider  a  multimode  seismic  wavetrain,  this  directional 
scattering  will  work  to  reduce  the  reinforcement  of  modal  amplitudes  by  cross¬ 
coupling.  An  average  of  the  estimated  scattering  attenuation  across  a  number  of 
modes  will  therefore  give  a  measure  of  the  decay  rate  of  the  wavetrain  with  distance 
which  can  be  compared  with  observations  on  attenuation.  It  should  also  be  recalled 
that  the  derivation  of  the  scattering  attenuation  depends  on  propagation  occurring 
over  distances  which  are  long  compared  with  the  correlation  length  of  the 
heterogeneity. 
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Estimates  of  Scattering  Attenuation 

The  computational  procedure  which  we  have  just  derived  enables  us  to  generate 
a  good  estimate  of  the  statistical  loss  factor  for  an  individual  mode  as  a  combination 
of  anelastic  and  scattering  contributions.  For  direct  comparison  between  different 
modes  we  can  use  the  rate  of  energy  decay  with  distance.  We  will  here  apply  our 
results  to  the  regional  phases  Lg,  Sn  travelling  through  heterogeneous  earth  models. 

From  equation  (25),  calcubtion  of  the  scattering  loss  factor  requires  us  to  specify 
the  statistical  nature  of  the  horizontal  heterogeneity  spectrum  in  the  crust  and  upper 
mantle  as  well  as  estimate  the  ensemble  averages  of  the  intermode  coupling  terms  , 

which  depend  on  the  modal  eigenfunctions  for  the  reference  model.  In  figure  1 
we  show  our  reference  model  (ARANDA)  together  with  the  vertical  component  of 
displacement  for  Rayleigh  modes  at  IHz.  The  modes  clearly  divide  into  two  classes; 
firstly,  those  modes  for  which  the  displacement  is  largely  restricted  to  the  crust  which 
will  constitute  the  Lg  wave  train,  and  secondly,  modes  with  little  crustal  displacement 
but  significant  energy  transport  in  the  mantle  which  represent  the  Sn  phase.  At  1  Hz 
the  transition  occurs  at  mode  11.  The  fundamental  and  first  higher  mode  are  confined 
to  the  sedimentary  overburden  and  hardly  interact  with  the  other  modes. 

The  ensemble  averages  require  an  average  of  depth  integrals  over  a 

combination  of  heterogeneity  and  eigenfunction  tenns.  We  have  estimated  these 
averages  by  a  Monte  Carlo  simulation:  random  velocity  perturbations  (up  to  1  per 
cent)  were  applied  to  the  velocities  at  the  top  and  bottom  of  each  of  the  layers  in  the 
reference  model  ARANDA  and  linear  interpolation  was  used  to  calculate  intermediate 
velocity  values.  Densitites  were  varied  in  proportion  to  the  velocities  and  the  depths 
of  the  interfaces  were  not  varied.  The  coupling  integrals  (3)  were  then  evaluated  for 
both  Love  and  Rayleigh  waves  for  300  different  simulations  of  the  vertical 
heterogeneity  structure  and  averaged  to  give  estimates  of  is  calcubted 

for  an  integral  over  the  crust  i.e.  down  to  30  km  and  over  the  mantle  structure 
which  is  extended  in  the  reference  model  to  200  km  depth.  This  procedure  gives  quite 
stable  estimates  for  and 

The  model  adopted  for  the  autocorrelation  of  the  horizontal  variation  of  the 
heterogeneity  was  exponential 

R{u)  =  /i’exp(-l»  I/D), 

for  a  horizontal  correlation  length  D.  This  form  has  the  advantage  that  it  can  allow 
the  existence  of  discontinuities  in  velocity  gradients.  The  corresponding  Fourier 
transfonn  appearing  in  the  loss  factor  terms  is 
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<\F  {k)\^  >  =  2h^/Dlkh{-[lD9]. 


We  have  scaled  the  amplitude  factor  h  to  be  unity  for  the  ±1  per  cent  heterogeneity 
assumed  in  the  calcubtion  of  the  vertical  averages  . 

The  dependence  of  the  loss  factor  on  heterogeneity  will  be  quadratic  over  the 
ranges  for  which  the  present  theory  is  valid  (less  than  ±2  per  cent).  At  higher  levels, 
local  multiple  scattering  will  become  important  and  the  rate  of  increase  of  the  loss 
factor  will  drop  below  quadratic. 

In  figure  2,  we  show  the  behaviour  of  the  scattering  loss  factor  for  three  Rayleigh 
and  Love  modes  at  1  Hz,  with  ±1  per  cent  heterogeneity,  as  a  function  of  the 
correlation  length  for  the  crust.  The  mantle  correlation  length  was  fixed  at  50  km. 
The  modes  were  chosen  to  represent  different  aspects  of  the  wavefield.  At  1  Hz,  mode 
4  is  sensitive  to  the  upper  and  middle  crust  and  forms  part  of  the  onset  of  Lg  with  a 
group  velocity  of  .3.45  km/s.  Mode  9  has  greater  sensitivity  to  lower  crustal  properties 
and  represents  waves  travelling  near  critical  incidence  on  the  crust-mantle  interface;  at 
1  Hz  it  is  close  to  an  Airy  Phase  with  a  group  velocity  of  3.2  km/s  and  so  represents 
the  coda  of  Lg.  Mode  14  has  a  group  velocity  around  4.45  km/s  and  forms  part  of  the 
Sn  phase. 

The  solid  triangles  in  figure  2  indicate  the  calcubtions  for  two-dimensional 
structure  for  Rayleigh  waves  and  the  open  diamonds  indicate  the  corresponding 
values  for  Love  waves.  The  behaviour  is  generally  similar  although  there  is  more 
difference  for  the  Sn  type  mode  (14).  For  this  slowness  the  Rayleigh  w’ave  couples  P 
waves  propagating  in  the  near  surface  with  Sn  type  behaviour  for  S;  the  influence  of 
the  P  wave  velocity  heterogeneity  raise  the  loss  factor  for  the  mode. 

For  both  of  the  crustal  mode  we  see  a  tendency  for  the  scattering  loss  factor  to 
peak  for  correlation  distances  about  3  times  the  horizontal  wavelength  of  the  waves. 
For  small  scale  heterogeneity  there  will  be  considerable  scattering  but  loss  and  gain  by 
intermodal  interactions  will  tend  to  balance.  For  larger  scale  heterogeneity,  scattering 
w’ill  be  more  infrequent  and  so  loss  will  be  reduced.  The  detailed  behaviour  depends 
on  the  mode  and  also  the  character  pre.scribed  for  the  vertical  heterogeneity. 

The  loss  factor  for  the  mantle  mode  (14)  is  only  less  sensitive  to  the  character  of 
the  crustal  heterogeneity.  There  is  considerable  mixing  between  mantle  modes  giving 
a  significant  scattering  component.  However,  the  intrinsic  attenuation  is  likely  to  be 
low  so  that  the  overall  loss  factor  will  normally  be  less  than  for  crustal  propagation. 

In  order  to  compare  a  number  of  different  modes  it  is  preferable  to  work  with 
the  rate  of  energy  loss  with  distance  h„  rather  than  the  scattering  loss  factors.  In 
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figure  3  we  therefore  display  b„  as  a  function  of  mode  number  at  IHz,  for  both  Love 
and  Rayleigh  waves  and  three  different  choices  of  crustal  correlation  length  D^.  Mode 
1  is  confined  to  the  sediments  and  has  little  interaction  with  other  modes,  so  that  its 
scattering  attenuation  is  small.  Attenuation  for  the  other  modes  is  of  the  same  order 
of  magnitude.  As  the  horizontal  correlation  length  shortens  to  less  than  10  km,  a 
reasonable  value  for  crustal  variation  (Wu  &  Aki  1988),  the  decay  rates  of  the  various 
crustal  modes  tend  to  equalise  which  will  give  rise  to  a  consistent  decay  rate  for  the 
Lg  wavetrain. 

The  frequency  dependence  of  the  toss  factor  estimates  is  explicit  through  the 
dependence  on  co,  k  in  the  coupling  coefficients  and  correlation  spectra,  but  also 
implicit  through  the  shape  of  the  eigenfunctions.  For  a  crustally  guided  mode  it  is 
difficult  to  compare  the  Q”’  estimates  at  different  frequencies  because  the  character  of 
the  mode  changes  However,  we  can  compare  modes  with  similar  propagation 
characteristics  at  different  frequencies  and  we  find  that  our  statistical  estimates  of  sQ"’ 
increase  with  frequency.  The  increase  is  not  a  simple  power  law  and  these  is  an 
indication  of  saturation  at  higher  frequencies  (around  2  Hz). 

As  pointed  out  by  Hudson  (1982)  it  is  often  difficult  to  find  a  direct  relation 
between  the  results  of  stochastic  calculations  and  observable  features  of  the  seismic 
wavefield.  Fortunately,  for  guided  wave  attenuation  we  are  able  to  find  a  very  close 
correspondence  between  our  estimate  of  energy  loss  with  distance  and  observations  of 
regional  phase  attenuation  over  a  variety  of  paths  within  a  single  region  (see  eg. 
Nuttli  1980,  Herrmann  1980).  Once  the  effects  of  dispersion  and  geometrical 
spreading  are  removed  from  the  observations,  the  attenuation  with  distance  is  isolated 
and  will  include  both  intrinsic  and  scattering  attenuation. 

For  a  large  number  of  different  patlrs  in  Iran  at  1  Hz,  Nuttli  (1980)  found  a  broad 
spread  of  amplitude  decay  factors  for  both  Lg  and  Sn  with  an  average  value  around 
0.0045  km“'.  To  convert  to  energy  loss  coefficients  we  must  mutiply  by  a  factor  of 
two.  Comparison  with  figure  3  shows  that  ±1  per  cent  heterogeneity  would  account 
for  about  10  per  cent  of  the  observed  attenuation  and  this  factor  would  rise  to  around 
40  per  cent  for  ±2  per  cent  heterogeneity.  The  intrinsic  attenuation  is  likely  to  be  high 
in  Iran,  but  the  apparent  Q  of  around  2(X)  may  well  have  a  significant  scattering 
component. 

For  the  low  attenuation  zone  in  the  Eastern  United  States,  Herrmann  (1980)  has 
dedueexj  an  amplitude  decay  rate  of  0.0009  km  ',  with  an  apparent  Q  of  1500,  for 
many  paths  to  the  station  BLA  (Black.sburg,  Virginia).  This  would  be  equivalent  to  an 
energy  decay  rate  of  0.(X)18  km  Nearly  a  third  of  this  loss  could  be  accounted  for 
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by  scaHering  due  to  ±1  per  cent  heterogeneity,  within  a  medium  with  low  intriiisic 
loss  (Q,  22'^0). 

These  results  show  that  our  estimate  of  the  scattering  loss  for  guided  waves  fit 
well  with  observations  of  guided  wave  attenuation.  If,  then,  we  have  a  measure  of 
the  heterogeneity  structure  and  intrinsic  attenuation  for  a  region  we  can  use  our 
statistical  approach  to  estimate  and  average  effective  loss  factor  which  can  be  used  as 
a  reference  against  which  to  compare  observations  on  different  paths. 

The  seismic  velocity  distributions  derived  from  tomographic  inversion  (e.g. 
Spakman  1989)  often  show  higher  levels  of  heterogeneity  than  have  been  assumed  in 
the  stochastic  calculations  (up  to  ±5  per  cent).  However  such  heterogeneity  has  a 
horizontal  scale  length  typically  greater  than  200  km  and  so  will  not  be  a  significant 
contributor  to  scattering  for  higher  frequency  waves  (1  Hz  and  above).  The  estimates 
of  scattering  attenuation  due  to  small  scale  variations  in  seismic  properties  will 
therefore  give  d  good  guide  to  the  behaviour  in  the  presence  of  variability  on  much 
longer  horizontal  scales. 
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figure  Captions: 

1.  ARANDA  reference  model  and  the  vertical  component  of  displacement  for  the 
first  18  Rayleigh  modes  at  1  Hz. 

2.  Variation  of  the  scattering  loss  factor  for  different  Rayleigh  and  Love 

modes  at  1  Hz  as  a  function  of  horizontal  correlation  length  in  the  crust. 
Rayleigh  waves  are  indicated  by  solid  triangles  and  Love  waves  by  open 
diamonds. 

3.  Variation  of  the  energy  loss  with  distance  h„  with  mode  number  for  different 
horizontal  correlation  lengths  in  the  crust.  Rayleigh  waves  are  indicated  by  solid 
triangles  and  Love  waves  by  open  diamonds. 
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FIGURE  3 
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SUMMARY 

The  Lg  phase  is  one  of  the  most  prominent  features  of  many  regional 
seismograms  and  its  amplitude  is  known  to  be  sensitive  to  major  changes  in 
crustal  structure.  This  sensitivity  can  be  exploited  to  investigate  the  nature  of 
crustal  heterogeneity,  but  also  reduces  the  utility  of  the  analysis  of  Lg  phases  for 
characterizing  seismic  sources. 

The  method  of  ray  diagrams,  based  on  the  constructive  interference  of 
multiple  S  reflections,  gives  a  visual  measure  of  the  interaction  of  guided  wave 
trains  with  variations  in  crustal  structure.  The  technique  is  extended  to  three 
dimensional  heterogeneity  and  by  representing  the  S  ray  systems  in  plan  view 
through  their  basement  reflection  points,  a  semiquantitative  measure  of  the 
disruption  of  Lg  at  different  phase  velocities  can  be  made.  In  addition,  the  tilt  of 
the  ray  systems  can  be  used  as  an  indication  of  the  degree  of  interconversion 
between  Rayleigh  and  Love  wave  energy. 

Two  simple  crustal  models  are  examined  to  investigate  the  effects  of  oblique 
incidence  and  transverse  topographic  gradients  on  multiple  S  reflections  using  ray 
diagrams.  The  method  is  then  applied  to  several  real  earth  models  of  the  crust  in 
central  Asia  using  topographic  data  and  assuming  isostasy.  Interpretations  of  the 
ray  diagrams  for  paths  across  the  major  tectonic  features  compare  well  with  many 
of  the  observations  from  previous  studies  (Ruzaikan  et  al.,  1977;  Ni  and  Barazangi, 
1983),  and  suggest  that  the  gross  characteristics  of  Lg  propagation  can  be  attributed 
to  the  shape  of  the  crustal  wave  guide. 

INTRODUCTION 

The  Lg  phase  propagates  in  the  earth's  crust  and  dominates  seismograms  at 
regional  distances.  It  is  characterized  by  frequencies  between  1  and  5  Hz  and  is 
often  of  extended  duration.  Although  the  wavetrain  does  not  exhibit  a  clear  onset, 
it  builds  to  an  a  amplitude  maximum  at  a  group  velocity  close  to  3.5  km/s  and 
may,  in  regions  of  sedimentary  cover,  carry  significant  energy  to  group  velocities  of 
2.8  km/s.  Recently  considerable  attention  has  focussed  on  the  signature  of  the  Lg 
wavetrain  in  deducing  crustal  morphology  and  as  an  aid  to  nuclear  discrimination; 
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in  both  instances  an  understanding  of  the  manner  in  which  Lg  interacts  with 
zones  of  crustal  heterogeneity  is  a  primary  objective. 

The  sensitivity  of  Lg  to  major  variation  in  crustal  structure  has  been  widely 
noted;  early  studies  revealed  that  as  little  as  100  km  of  intervening  oceanic  crust  is 
sufficient  to  eliminate  Lg  propagation  across  ocean  basins  (Ewing  et  al.,  1957).  The 
character  of  Lg  propagation  across  certain  continental  areas  where  the  crust  is 
known  to  be  complex,  most  notably  central  Asia,  has  been  used  to  constrain 
structure  and  provide  insight  into  the  operative  tectonic  processes  (Ruzaikan  et  al., 
1977;  Kadinsky-Cade  et  al.,  1981;  Ni  and  Barazangi,  1983).  These  studies  have 
involved  the  interpretation  of  short  period  seismograms  to  qualitatively 
characterize  broad  geographical  zones  in  terms  of  their  amenability  to  Lg 
propagation.  More  quantitative  analysis  is  possible  with  increased  raypath 
coverage;  Kennett  et  al.  (1985)  inverted  a  comprehensive  set  of  Lg  attenuation  data 
for  paths  in  NW  Europe  to  delineate  crustal  heterogeneity  in  the  North  Sea  region. 

Two  different  descriptions,  which  are  equivalent  in  the  case  of  a  stratified 
medium,  can  be  made  for  the  physical  nature  of  Lg  in  a  heterogeneous  crust:  the 
phase  can  be  considered  as  a  sum  of  higher  mode  surface  waves  whose  energy  is 
mostly  confined  to  the  crust,  or  alternatively,  as  the  result  of  constructive 
interference  between  S-waves  multiply  reflected  between  the  free  surface  and  the 
crust-mantle  boundary.  For  heterogeneous  zones,  two  methods  have  evolved  to 
model  the  propagation  of  Lg  based  on  these  different  views  of  the  propagation 
process. 

The  coupled  mode  technique  (Kennett,  1984)  is  based  on  the  modal 
interpretation  and  solves  the  boundary  value  problem  exactly;  however  in  zones 
of  exaggerated  heterogeneity  and  variable  layer  thickness  it  can  become 
computationally  intractable.  The  method  of  ray  diagrams  is  a  simpler  means  of 
assessing  the  propagation  of  Lg  in  a  semi-quantitatitive  manner  and  is  not  limited 
by  the  same  circumstances.  It  involves  tracing  a  set  of  rays  at  fixed  slowness 
through  a  crustal  layer  of  variable  thickness,  and  the  disruption  in  the  pattern  of 
the  ray  system  (or  equivalently  the  reflection  points)  is  used  as  a  measure  of  modal 
coupling  and  scattering.  Kennett  (1986)  compared  the  visually  descriptive  ray 
diagrams  with  coupled  mode  solutions  for  simple  2D  earth  models  and  established 
the  validity  of  the  ray  based  technique. 

Our  objective  in  this  study  is  to  extend  Kennett's  (1986)  approach  to  3D 
heterogeneity.  After  describing  a  basis  for  the  interpretation  of  ray  diagrams  we 
examine  two  simple  earth  models  and  investigate  the  effects  of  oblique  incidence 
on  Lg  propagation.  We  then  shift  our  focus  to  a  central  Asia  to  compare  our 
predictions  of  Lg  propagation  with  observations  from  previous  studies. 
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THE  RAY  DIAGRAM  METHOD 

Kennett  (1986)  provided  a  basis  for  the  interpretation  of  ray  diagrams  by 
comparing  them  with  coupled  mode  solutions  for  a  suite  of  2D  earth  models  and 
noting  some  simple  and  consistent  relationships.  The  coupled  mode  technique 
describes  Lg  as  a  sum  of  modal  eigenfunctions,  associated  with  a  laterally 
homogeneous  reference  model,  and  weighted  by  modal  coefficients  which  vary 
w'ith  position.  In  2D  these  coefficients  satisfy  a  system  of  non-linear  differential 
equations  which  can  be  solved  exactly  in  principle.  In  practice,  however,  the 
method  becomes  computationally  intractable  in  zones  of  exaggerated  heterogeneity 
(greater  than  5%  velocity  perturbation  or  2  km  in  boundaiy  perturbation  for  1  Hz 
waves)  as  consideration  of  intermode  coupling  among  an  increasing  number  of 
modes  becomes  necessary.  Maupin  (1988)  has  made  use  of  coupled  local  modes  to 
give  a  more  direct  treatment  of  variable  boundaries;  however  changes  in  position 
of  interfaces  require  extensive  recomputation  of  modal  eigenfunctions.  Ray 
diagrams,  in  contrast,  are  most  useful  where  the  shape  of  the  crustal  waveguide 
varies  and  become  more  difficult  to  interpret  if  a  non-uniform  velocity  structure  is 
imposed.  The  analysis  of  Kennett  (1986)  included  models  amenable  to  both 
methods,  and  demonstrated  the  close  correspondence  between  the  coupled  mode 
and  ray  diagram  results. 

The  correspondence  between  the  two  methods  is  most  easily  explained  by 
noting  that  for  a  given  frequency  the  Lg  wave  train  will  comprise  a  finite  number 
of  modes.  Each  mode  is  characterized  by  a  specific  phase  velocity  which  can  be 
associated  with  an  S-wave  ray  system  at  a  particular  angle  to  the  vertical.  This 
suggests  that  the  geometrical  regularity  of  a  system  of  rays  launched  in  some 
systematic  fashion  may  be  used  as  an  indication  of  the  coherence  and  sustained 
amplitude  of  a  given  mode  in  addition  to  a  measure  of  the  S-wave  constructive 
interference  condition. 

At  fixed  frequency,  the  angle  of  S-wave  propagation  to  the  vertical  increases  as 
the  phase  velocity  decreases.  For  a  given  angle  to  the  vertical,  determined  by  the 
choice  of  incident  mode,  the  behaviour  after  interaction  with  heterogeneity  can  be 
conveniently  characterized  by  the  spread  in  propagation  angles.  Increased  angles  to 
the  vertical  correspond  to  conversion  to  lower  order  modes  with  the  change  in 
angle  proportional  to  the  spread  in  intermode  coupling.  Similarly  steeper 
propagation  corresponds  to  conversion  to  higher  order  modes.  In  general  waves 
incident  with  higher  phase  velocities  tend  to  show'  a  larger  variation  in  ray  angle 
after  transmission  through  the  heterogeneity  due  to  increased  multiple  reflections 
in  the  complex  zone.  In  those  cases  where  caustics  develop  in  the  ray  diagrams  the 
equivalent  coupled  mode  solutions  show  considerable  variation  in  the  amplitudes 
of  the  modal  coefficients  with  position. 
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The  crustal  models  presented  in  the  following  analysis  exhibit  3D 
heterogeneity.  It  is  no  longer  possible  to  construct  ray  patterns  from  a  system  of 
plane  waves  in  cross  section  as  in  the  2D  case.  Rather,  we  have  resorted  to  plan 
views  representing  surface/Moho  topography  via  contours  and  consider  a  single 
point  source  emitting  rays  at  a  given  angle  to  the  vertical  (modal  slowness) 
through  a  range  of  azimuths.  Regularity  in  plotted  surface/Moho  reflection  points 
and  the  horizontal  projection  of  raypaths  provide  a  measure  of  the  constructive 
interference  condition  in  a  manner  analogous  to  the  pattern  of  2D  ray  systems 
considered  by  Kennett  (1986).  An  added  feature  of  the  3D  analysis  is  that  we  are 
able  to  examine  the  way  in  which  the  wavefronts  tilt  as  they  encounter  features 
oblique  to  the  line  of  propagation.  This  'ray  tilt'  will  lead  to  a  rotation  of  the 
polarization  of  the  S  waves  between  the  vertical  and  horizontal  planes.  As  a 
result  energy  having  an  explosive  source  lying  purely  in  a  vertical  plane  (i.e. 
represented  solely  by  higher  mode  Rayleigh  waves),  may  after  encounter  with  the 
heterogeneity  have  some  component  in  the  horizontal  plane,  which  would  be 
interpreted  as  conversion  to  higher  mode  Love  waves.  Zones  of  the  model  where 
the  ray  patterns  exhibit  a  consistent  change  in  tilt  can  be  expected  to  produce 
significant  conversion. 

CRUSTAL  MODELS 

We  have  investigated  Lg  propagation  in  two  simple  models  of  crustal  structure 
and  in  reconstructions  of  the  crustal  structure  in  central  Asia  based  on  world 
topography  data,  using  ray  diagrams.  All  of  these  models  are  characterized  by 
simple  variation  in  surface  and  basement  (Moho)  topography;  physical  properties 
remain  constant  throughout  the  crust. 

Although  a  more  sophisticated  analysis  incorporating  crustal  stratification 
might  be  implemented  with  slight  modification,  this  has  not  been  done  for  several 
reasons.  Ray  diagrams  are  at  best  semiquantitative  and  are  meant  to  provide  a 
simple  means  of  assessing  Lg  propagation  through  crustal  structure.  Moreover, 
although  the  shape  of  the  crustal  waveguide  is  the  dominant  factor  influencing  Lg, 
other  factors  such  as  small  scale  heterogeneity  in  crustal  velocities  cannot  be 
entirely  ignored.  In  addition,  when  we  construct  a  crustal  model  from  just  the 
surface  topography  we  are  forced  to  make  a  number  of  assumptions  which  may  not 
he  well  founded,  and  so  a  very  complex  model  is  not  warranted. 

Shear  velocities  of  3.5  km/s  and  4.6  km/s  were  chosen  to  characterize  the  crust 
and  upper  mantle.  For  a  given  phase  velocity  rays  were  traced  away  from  a  point 
source  with  multiple  reflections  at  the  surface  and  the  Moho.  The  character  of  the 
propagation  patterns  are  indicated  by  plotting  a  plan  section  of  the  rays  and 
marking  the  reflection  points  at  the  Moho.  The  magnitude  of  changes  in  ray  tilt 


24 


are  illustrated  on  separate  ray  diagrams  by  'tics'  located  perpendicular  to  the  ray  at 
surface  and  Moho  reflection  points.  This  display  enables  the  regular  character  of 
the  Lg  wavefront  in  the  homogeneous  regions  to  be  compared  with  the  disrupted 
pattern  after  passage  through  the  heterogeneity. 

SIMPLE  CRUSTAL  STRUCTURES 

In  this  section  we  consider  two  structures  for  which  the  surface  and  Moho 
topography  are  represented  as  piecewise  smooth  surfaces  and  the  zones  of 
variation  are  of  limited  spatial  extent.  The  ray  patterns  resulting  from  these 
simple  configurations  are  intended  to  demonstrate  the  effects  of  simple  changes  in 
topography  and  aid  in  the  interpretation  of  more  complicated  ray  diagrams  for 
central  Asia  presented  in  the  following  section. 

Crustal  Thickening 

The  first  case  is  of  a  linear  mountain  chain  and  a  section  across  the  model 
along  the  profile  A-B  is  shown.  Ray  diagrams  were  examined  for  the  full  range  of 
phase  velocities  associated  with  Lg  propagation  (3.5  to  4.6  km/s)  and  several 
general  features  were  observed  to  characterize  the  wavefield  emergent  from  the 
structur..  In  figure  1  we  show  the  ray  diagrams  for  phase  velocities  of  3.7  km/s 
and  4.3  km/s  which  illustrate  the  major  features  of  the  wavefield.  Perpendicular 
to  the  axis  of  the  range  we  note  a  central  corridor  marked  by  a  continuous  and 
regular  arrangement  of  bottom  reflection  points  and  little  variation  in  ray  density 
with  azimuth  (even  though  the  angle  to  the  vertical  is  increased).  This  corridor  is 
flanked  by  zones  of  alternating  high  and  low  ray  density  as  the  angle  between  the 
axis  and  the  wavefront  increases. 

For  particular  combinations  of  source  location  and  phase  velocity,  some  zones 
may  exhibit  regularity  in  both  ray  density  and  reflection  pattern,  and,  therefore 
represent  windows  through  which  a  significant  fraction  of  the  initial  Lg  wave 
energy  can  pass.  As  ray  incidence  becomes  increasingly  oblique  to  the  mountain 
range  axis,  rays  become  trapped  within  the  range  and  few  rays  penetrate  the 
structure  to  the  crust  outside.  Thus  we  might  expect  linear  zones  of  crustal 
thickening  to  act  as  lateral  waveguides,  especially  for  earthquakes  in  close 
proximity. 

Rays  which  impinge  on  the  surface  at  angles  to  the  surface  normal  less  than 
the  critical  angle  defined  by  Snell's  law  will  undergo  transmission  into  the  upper 
mantle  and  are  distinguished  by  filled  diamonds  marking  the  Moho  reflection 
points  .  Kennett  (1986)  notes  that  this  situation  represents  a  coupling  to  modes 
constituting  the  Sn  phase.  As  observed  in  the  equivalent  2D  model  (Kennett, 
1986),  refraction  after  emergence  from  the  mountain  range  occurs  most 
prominently  near  the  maximum  phase  velocity  for  Lg  (above  4.4  km/s).  The 
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distribution  of  these  emergent  refracted  rays  is  dependent  on  source  location  and 
displays  little  azimuthal  preference.  In  figure  2  we  show  the  polarization  of  the  S 
rays  for  the  same  phase  velocities  as  figurel  by  means  of  'tics'  which  show  the 
projection  onto  a  horizontal  plane  of  waves  which  left  the  source  with  energy 
confined  to  the  vertical  plane.  Figure  2  reveals  that  changes  in  ray  tilt  develop  in 
the  zones  flanking  the  central  corridor  where  topographic  gradients  become 
increasingly  transverse  to  the  direction  of  propagation.  Thus  we  expect  some 
transfer  of  energy  between  vertical  and  horizontal  planes  which  would  be 
interpreted  as  a  conversion  between  Love  and  Rayleigh  waves.  This  feature 
introduces  further  complexity  to  what  is  already  seen  to  be  a  very  complicated 
wavefield. 

Crustal  Thinning 

The  second  case  represents  a  region  of  crustal  thinning  simulating  a  zone  of 
localized  crustal  extension  and  is  designed  to  illustrate  the  effect  of  local  transverse 
gradients  in  topography.  In  cross  section  it  is  characterized  by  flanks  rising  from  30 
km  to  20  km  depth  to  a  plateau  at  the  Moho;  its  configuration  in  the  horizontal 
plane  is  elliptical.  Kennett  (1986)  noted  in  an  analogous  2D  model  that  constriction 
of  the  crustal  waveguide  had  a  more  severe  effect  on  transmission  of  Lg  than 
increase  in  thickness.  This  is  also  apparent  in  the  3D  model  (see  figure  3)  where 
we  again  display  ray  patterns  for  phase  velocities  of  3.7  and  4.3  km/s.  The  effect  of 
the  crustal  thinning  is  to  steer  rays  away  from  the  angular  window  defined  by  the 
lateral  boundaries  of  the  structure  and  the  source  position.  This  effect,  due  to  the 
presence  of  strong  gradients  in  topography  oblique  to  the  ray  path,  gives  rise  to 
caustics  in  the  ray  pattern,  and  is  particularly  noticeable  at  higher  phase  velocities 
as  in  figure  3c.  The  constriction  of  the  waveguide  can  also  lead  to  waves  travelling 
across  the  top  of  the  structure  being  refracted  into  the  mantle,  and  hence  the  low 
density  of  emergent  rays  is  complemented  by  an  additional  energy  loss  through 
conversion  to  Sn. 

In  figure  4,  we  show  the  ray  tilt  representation  of  the  S  wave  polarisation  for 
this  crustal  thinning  model.  The  three  dimensional  character  of  the  crustal 
heterogeneity  leads  to  polarization  changes  for  a  larger  proportion  of  the  rays  than 
for  the  linear  mountain  chain  of  figure  2.  At  the  higher  phase  velocities  the 
polarization  rotation  is  particularly  strong  at  the  lower  margin  of  the  heterogeneity 
where  the  gradients  in  Moho  topography  are  nearly  transverse  to  the  ray  path. 

CENTRAL  ASIA 

The  propagation  of  Lg  is  of  considerable  interest  in  determining  the  nature  of 
the  earth's  crust  in  central  Asia.  Extreme  variations  in  topographic  relief  over  the 
region  shown  in  figure  .5  testify  to  a  complex  and  heterogeneous  waveguide.  To  the 
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north  of  the  Indian  shield  (average  elevation  200  m),  the  Himalayas  rise  above 
6000  m  and  form  the  southern  perimeter  of  the  Tibetan  Plateau,  a  remarkable 
physiographic  feature  covering  some  2.5  x  10^  km^  and  characterized  by  an  average 
elevation  of  about  5000  m.  The  Tarim  Basin  separates  the  Tibetan  Plateau  to  the 
northwest  from  a  second  major  mountain  range,  the  Tien  Shan. 

The  pattern  of  Lg  propagation  varies  significantly  over  the  region  (see  Figure  5; 
Ruzaikan  et  al.,  1977;  Ni  and  Barazangi,  1983).  Specifically,  Lg  propagates  very 
efficiently  across  the  more  tectonically  stable  areas:  the  Indian  Shield,  Tarim  Basin 
and  Eurasian  Platform.  The  phase  is  also  present  for  paths  along  the  Himalayas 
and  Tien  Shan  (although  less  clear);  however  it  is  generally  weak  or  absent  for 
paths  crossing  the  Tibetan  Plateau.  These  observations  bear  important 
implications  for  the  interpretation  of  the  tectonic  evolution  of  the  region. 
Additional  interest  in  the  nature  of  Lg  propagation  across  central  Asia  is  prompted 
by  the  presence  of  both  Soviet  and  Chinese  nuclear  test  sites.  The  prominence  of  Lg 
waves  on  regional  seismograms  makes  the  phase  valuable  for  the  detection  of 
small  events,  but  the  sensitivity  of  Lg  to  the  variation  in  crustal  structure  across 
the  region  reduces  its  utility  for  discrimination  and  a  better  knowledge  of  the 
effects  of  heterogeneity  on  propagation  is  necessary. 

Surface  topography  for  the  following  models  of  the  crustal  structure  was 
constructed  by  smoothing  (25  point  average)  digital  elevation  data  supplied  at  5' 
intervals  from  a  world  topography  data  base  (ETOP05  compiled  by  the  National 
Geophysical  Data  Center,  Boulder,  Colorado).  The  depth  to  the  Moho  was 
determined  by  further  smoothing  (81  point  average)  of  the  surface  topography  and 
assuming  90%  isostatic  compensation  occurs  within  the  crust.  This  assumption, 
although  simplistic,  is  thought  to  be  reasonable  within  the  limitations  of  the  ray 
method  for  a  uniform  crust. 

The  choice  of  a  suitable  phase  velocity  in  representing  Lg  propagation  in  the 
real  earth  by  ray  diagrams  deserves  some  discussion.  In  typical  crustal  situations 
the  higher  modes  constituting  Lg  propagate  at  a  range  of  phase  velocities.  Within 
this  range  certain  modes  within  loosely  defined  phase  velocity  windows  will 
dominate  at  different  group  velocities.  The  onset  of  Lg  on  short  period 
seismograms  is  determined  by  modes  at  phase  velocities  of  3.4  to  3.6  km/s  (i.e. 
crustal  shear  velocities).  The  amplitude  maximum  which  follows  is  usually 
characterized  by  slightly  greater  phase  velocities  between  3.7  and  4.0  km/s  whereas 
typical  values  for  the  late-arriving  Airy  phases  reach  4.3  to  4.5  km/s.  These  latter 
phases  play  an  important  role  in  determining  the  character  of  regional 
seismograms  especially  in  zones  of  heterogeneity  where  conversion  to  Sn  is 
•significant.  An  additional  factor  to  consider  is  the  effect  of  velocity  variation 
within  the  crust.  Although  the  shape  of  the  crustal  wave  guide  is  a  major 
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influence  on  Lg  propagation,  a  positive  crustal  velocity  gradient  in  depth  will  have 
some  effect  on  the  configuration  of  the  ray  systems,  flattening  them  and  extending 
the  distances  between  reflection  points.  In  our  simple  single  layer  crustal  model 
this  is  equivalent  to  a  decrease  in  phase  velocity.  By  examining  the  horizontal 
distance  between  reflections  for  single  layer  and  multilayered  crustal  models  based 
on  equivalent  total  traveltime,  we  note  that  the  actual  magnitude  of  the  decrease  is 
only  significant  at  very  low  phase  velocities  (i.e.  within  a  few  percent  of  the  crustal 
shear  velocity).  Therefore  in  modelling  effects  on  the  early  arrivals  it  may  be 
appropriate  to  compensate  by  reducing  the  phase  velocity  parameter  for  crustal 
models  of  uniform  velocity.  As  the  overall  character  of  the  Lg  wavetrain  is 
dependent  on  modal  contributions  at  a  full  range  of  phase  velocities  we  have 
attempted  to  select  model  values  which  both  illustrate  the  general  disruption  in  Lg 
propagation  due  to  changes  in  the  shape  of  the  crustal  wave  guide,  and  are 
pertinent  to  observations  made  in  previous  studies. 

Although  the  major  topographic  features  in  the  area  exist  at  a  larger  lateral 
scale  than  those  in  the  two  simple  structural  models,  the  ray  patterns  can 
nonetheless  be  quite  sensitive  to  source  position.  Attempts  have  been  made  in  the 
following  analysis  to  restrict  interpretation  to  those  features  which  remain  stable 
over  small  changes  in  source  location. 

Semipalitinsk 

For  our  first  study,  we  examine  Lg  propagation  in  the  vicinity  of  the  Soviet 
nuclear  test  sites  located  in  eastern  Kazakhstan,  approximately  400  km  northeast  of 
Lake  Balkhash,  near  Semipalitinsk  Seismic  stations  at  Novosibirsk  (NSB)  and 
Talgar  (TLG)  are  shown  as  triangles  to  the  north  and  southeast  of  the  source  at 
Semipalitinsk  in  figure  6.  The  rather  moderate  and  consistent  topography  over 
much  of  this  region  permits  an  a.ssessment  of  our  modelling  procedure  in  what  is 
presumably  a  relatively  uniform  portion  of  the  crustal  waveguide. 

The  ray  diagrams  showm  in  figure  6  are  generated  at  a  phase  velocity  of  4.3 
km/s  since  anomalous  propagation  patterns  are  most  pronounced  at  velocities 
above  4.0  km/s.  The  character  of  the  patterns  can  be  discussed  in  tenns  of  three 
broad  geographic  zones.  To  the  north,  elevation  decreases  very  gradually, 
consequently  the  ray  patterns  are  very  regular  and  experience  a  slight  increase  in 
phase  velocity  with  distance  from  the  source.  For  initial  phase  velocities  above 
4.45  km/s,  this  results  in  energy  loss  through  conversion  to  Sn,  however  Lg 
transmission  through  this  zone  for  lower  phase  velocities  is  likely  to  be  extremely 
efficient.  It  is  worth  noting  that  Lg  recorded  at  NORSAR  would  propagate  through 
this  region.  The  general  topography  in  Kazakhstan  to  the  west  of  the  test  sites 
remains  relatively  constant  although  there  are  minor  local  fluctuations.  Ray 
patterns  in  this  zone  are  slightly  less  coherent  but  still  suggest  efficient  Lg 
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transmission.  The  Tien  Shan  and  Altai  chains,  mountain  ranges  rising  over  5000 
m,  border  the  region  to  the  .south  and  east.  Changes  in  relief  are  extreme  and  ray 
patterns  are  severely  distorted  especially  at  higher  phase  velocities.  The  character 
of  Lg  recorded  at  stations  along  these  paths  beyond  the  mountainous  regions  will 
undoubtedly  be  very  complex;  there  is  likely  to  be  considerable  variability  with 
position  and  significant  attenuation  due  to  conversion  to  Sn  and  scattering. 

In  the  corresponding  ray  till  diagram  (figure  6b),  as  expected,  rays  passing 
through  mountainous  areas,  where  elevation  gradients  transverse  to  direction  of 
propagation  are  likely  to  be  most  pronounced,  show  the  greatest  degree  of  tilt.  Also 
of  interest  and  perhaps  of  greater  significance  to  the  interpretation  of  source 
parameters  from  regional  seismograms  are  the  tilt  patterns  to  the  north  and  west. 
The  slight  but  relatively  uniform  north-south  gradient  in  elevation  over  these 
areas  causes  a  very  regular  and  obvious  tilt  pattern  in  eastern  Kazakhstan;  whereas 
north  of  Semipalitinsk,  where  rays  lie  subparallel  with  the  gradient,  virtually  no 
tilt  is  evident.  This  implies  that  Lg  recorded  at  stations  within  these  two  distinct 
windows  might  differ  markedly  in  the  relative  energy  distribution  among 
tran.sverse,  radial  and  vertical  components.  An  explosive  source  at  Semipalitinsk 
generating  shear  energy  primarily  on  vertical  and  radial  components  might 
nonetheless  give  rise  to  Lg  in  south  and  central  Europe  w'ith  a  dominant 
transverse  component,  while  Lg  recorded  further  north,  at  NORSAR  for  example, 
would  retain  much  of  its  Rayleigh-type  energy. 

Eastern  Tien  Shan 

We  now  turn  our  attention  to  a  more  complicated  region  comprising  the 
Tarim  Basin,  the  Altai  chain,  much  of  the  Tien  Shan  and  the  northern  fringe  of 
the  Tibetan  Plateau.  Ruzaikan  et  al.  (1977)  noted  that  earthquakes  in  the  eastern 
Tien  Shan  radiate  clear  Lg  phases  to  station  NSB  whereas  at  TLG  the  phase  is  often 
very  weak. 

Figures  7a  displays  the  reflection  patterns  at  phase  velocities  of  3.7  and  4.3 
km/s  for  a  source  in  this  area  of  interest.  The  low-er  velocity  diagram  reveals  a 
more  consistent  and  uniform  ray  pattern  at  NSB  than  TLG  even  though  the 
former  station  is  some  500  km  further  from  the  source.  This  observation  would 
tend  to  support  the  explanation  of  Ruzaikan  et  al.  that  the  difference  between 
signal  at  1  LG  and  NSB  results  from  different  path  lengths  within  the  complicated 
structure  of  the  Tien  Shan.  The  higher  phase  velocity  used  to  generate  the 
diagram  in  figure  7b  has  produced  a  more  severely  distorted  pattern  due  to  an 
increased  number  of  multiple  reflections.  Several  specific  points  are  worthy  of 
note.  Although  the  ray  density  and  patterns  are  still  somewhat  consistent  at  NSB, 
rays  in  the  vicinity  of  TLG  are  funnelled  into  two  arms  of  the  Tien  Shan  in  a 

manner  similar  to  that  de.scribed  for  the  simple  crustal  thickening  model.  Rays  in 
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the  northern  subsidiary  arm  rapidly  lose  energy  through  transmission  to  the 
mantle  whereas  those  in  the  southern  arm  remain  trapped.  This  funnelling 
behaviour  first  becomes  apparent  at  a  phase  velocity  of  3.9  km/s  and  suggests  that 
the  main  arm  may  act  as  a  lateral  wave  guide,  trapping  and  guiding  ^'Pergy  within 
a  zone  of  relatively  narrow  width  (100-200  km).  The  observations  of  Ruzaikan  et  al. 
would  tend  to  suggest  that  the  station  at  TLG  may  lie  just  outside  the  borders  of 
such  a  guide. 

Ni  and  Barazangi  (1983)  note  that  Lg  propagation  in  the  Himalayas  is  less 
efficient  across  strike  than  along  strike.  This  and  studies  in  other  mountainous 
regions  (e.g.  along  the  Andes  [Chinn  et  al.,  1980])  suggest  that  mountain  chains 
may  play  an  important  role  in  influencing  patterns  of  Lg  propagation  by  focussing 
and  guiding  energy  along  quasi-linear  zones  of  increased  crustal  thickness. 

Kunlun 

Lg  propagation  across  the  Tibetan  Plateau  has  been  the  focus  of  many  studies 
into  the  nature  of  the  crust  near  the  continental  collision  margin  Several  authors 
have  noted  the  absence  of  Lg  for  paths  crossing  the  plateau  and  the  peculiar  pattern 
of  propagation  for  sources  within  and  on  its  boundaries.  In  the  latter  case  there 
appear  to  be  well  defined  geographic  boundaries  to  effective  propagation,  notably 
the  36°  line  of  latitude  to  the  north  (Ruzaikan  et  al.,  1977)  and  the  litdus  Tsangpo 
suture  to  the  south  (Ni  and  Barazangi,  1983).  Lg  is  rarely  recorded  from  sources 
within  these  boundaries. 

Ray  diagrams  at  phase  velocities  of  3.7  and  4.3  km/s  for  a  source  just  above  the 

northern  boundary  in  the  Kunlun  are  shown  in  figure  8.  At  higher  phase 

velocities  (figure  8b),  a  large  proportion  of  rays  are  confined  to  the  plateau  by  lateral 

reflection  along  its  borders,  and  those  which  do  escape  have  undergone  so  many 

vertical  reflections  at  subcritical  angles  as  to  be  of  no  significance.  In  contrast, 

diagrams  at  lower  phase  velocities  maintain  regular  ray  patterns  through  most  of 

the  area.  Certainly,  patterns  are  fairly  regular  across  the  Tarim  Basin  where 

propagation  is  observed  to  be  reasonably  efficient  (Ruzaikan  el  al.,  1977).  Ray 

systems  emerging  from  the  Himalyas  are  also  fairly  coherent  hut  lose  energy 

through  transmission  to  the  mantle.  Our  as.sumptions  regarding  the  structure  of 

the  Moho  are  no  longer  valid  in  this  region  (see  Molnar,  1988),  where  the  crust  is 

quite  certainly  not  isostatically  compensated;  however  a  rapid  transition  in  crustal 

thickness  does  exist  and  should  complicate  Lg  propagation  in  a  similar  fashion. 

Efficient  propagation  across  the  Tarim  Basin  for  this  particular  source  location 

ceases  above  a  phase  velocity  of  approximately  4.0  km/s  where  the  effects  of  source 

position,  with  phase  velocity,  relative  to  the  crustal  transition  zone  (di.scussed  by 

Kennett  [1986])  become  manifest.  The  ray  diagram  approach  agrees  with 

suggestions  of  previous  authors  (Ruzaikan  et  al.,  1977)  that  the  absence  of  Lg  along 
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paths  crossing  the  Tibetan  Plateau  is  due  to  changes  in  the  crustal  waveguide, 
notably  constriction  along  its  margins.  Discriminating  between  natural  and 
artificial  sources  from  locations  within  Tibet  using  Lg  is  rendered  exceedingly 
difficult,  if  not  impossible,  by  this  propagation  effect. 

DISCUSSION 

We  have  seen  how  the  simple  device  of  ray  diagrams  can  be  extended  to  three 
dimensions  and  used  to  provide  insight  into  the  nature  of  the  interaction  of  the  Lg 
phase  with  heterogeneity.  The  interpretation  of  ray  diagrams  for  crude 
constructions  of  the  crust  in  central  Asia,  based  on  isostatic  compensation,  agree 
with  observations  from  previous  studies  and  indicate  that  much  of  the  general 
character  of  Lg  propagation  can  be  attributed  to  changes  in  the  shape  of  the  crustal 
waveguide.  Zones  of  regular  reflection  patterns,  such  as  those  in  Kazakhstan  for 
sources  at  Semipalitinsk,  characterize  a  waveguide  of  uniform  thickness  and 
generally  correspond  to  zones  of  efficient  Lg  propagation.  Rapid  decreases  in 
crustal  thickness,  as  for  example  the  margins  of  the  Tibetan  Plateau,  form 
boundaries  which  severely  restrict  Lg  propagation  especially  at  higher  phase 
velocities.  Where  such  changes  occur  over  narrow,  quasi-linear  zones  of  crustal 
thickening  (mountain  ranges  such  as  the  Tien  Shan,  or  the  Himalayas),  the 
possibility  exists  for  the  development  of  lateral  wave  guides.  The  three 
dimensional  ray  representation  also  affords  a  simple  means  of  assessing  the 
transfer  of  energy  between  horizontal  and  vertical  planes  through  examination  of 
the  change  in  ray  tilt  with  propagation.  Over  areas  with  transverse  gradients  in 
surface  and  Moho  topography  oblique  to  the  direction  of  propagation,  we  expect 
there  to  be  significant  energy  transfer  corresponding  to  an  interconversion  between 
Love  and  Rayleigh  modes.  This  feature  may  be  of  some  value  in  separating 
propagation  and  source  effects  in  the  problem  of  nuclear  discrimination  where 
energy  at  the  explosive  source  is  predominantly  F-5V. 

By  incorporating  detailed  information  on  gravity  anomalies  and  crustal 
structure  where  these  are  available  it  may  be  possible  to  further  improve  the 
accuracy  of  the  technique.  However,  caution  should  be  exercised  not  to  over 
extend  the  limits  of  its  usefulness;  internal  crustal  structure  undoubtedly  plays  a 
significant  and  complicated  role  in  the  propagation  of  Lg  as  is  made  evident  in  its 
effect  on  dispersion.  The  value  of  the  ray  approach  lies  in  assessing  the  main  effects 
of  crustal  shape  on  Lg. 
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FIGURE  CAPTIONS 


Figure  1 .  Reflection  patterns  for  a  linetir  zone  of  crustal  thickening,  simulating  a  mountain 
range.  A  vertical  cross  section  through  A-B  is  shown  in  a),  while  b)  and  c)  show  ray 
diagrams  generated  at  phase  velocities  of  3.7  and  4.3  km/s.  Moho  reflection  points  are 
plotted  as  open  squares;  filled  diamonds  indicate  leakage  out  of  the  crustal  waveguide. 

Figure  2.  Tilt  patterns  for  the  same  configuration  as  figure  I .  The  effects  of  transverse 
gradients  in  topography  are  evidenced  through  the  development  of  ray  tilt  at  angles  oblique 
to  the  axis  of  the  zone  of  thickening.  The  magnitude  of  a  tick  correspionding  to  a  90° 
change  in  tilt  is  indicated. 

Figure  3.  Reflection  patterns  impinging  on  an  elliptical  zone  of  crustal  thinning. 

Figure  4.  Tilt  patterns  for  the  same  configuration  as  figure  3. 

Figure  5.  Central  Asia,  a)  Topography  and  epicentral  data  from  Ruzaikan  et  al.  (1977) 
illustrating  characteristics  of  Lg  propagation  to  Talgar  (TLG)  are  plotted  as  circles  (Lg 
clear),  diamonds  (Lg  weak),  and  triangles  (Lg  absent),  b)  Reference  map  identifying  major 
topographic  features  in  a)  and  relative  postitions  of  maps  in  figures  6, 7  and  8.  Note  that 
there  is  some  systematic  distortion  in  this  and  the  following  figures  due  to  conversion  from 
spherical  to  rectangular  coordinates. 

Figure  6.  Reflection  (a)  and  tilt  (b)  patterns  for  a  source  at  the  Soviet  nuclear  test  site  near 
Semipalitinsk  at  a  phase  velocity  of  4.3  km/s.  Rays  in  the  northwestern  quadrant  represent 
energy  propagating  toward  western  Europe  and  show  very  regular  patterns  indicative  of 
efficient  Lg  propagation.  TTie  slight  but  regular,  nonh-south  gradient  in  topography  causes 
rays  over  the  western  half  of  the  map  in  (b)  to  exhibit  significant  tilt,  indicating  a  transfer  of 
energy  between  horizontal  and  vertical  planes.  The  magnitude  of  a  tick  corresponding  to  a 
90°  change  in  tilt  is  indicated. 

Figure  7.  Reflection  patterns  for  a  source  in  the  eastern  Tien  Shan  at  phase  velocities  of 
(a)  3.7  km/s  and  (b)  4.3  km/s.  Note  distortion  in  reflection  pattern  eastwards  toward 
Talgar  in  (a).  Energy  is  guided  into  two  arms  of  the  Tien  Shan  near,Talgar  in  (b),  a  feature 
which  first  becomes  apparent  at  a  phase  velocity  of  3.9  km/s. 

Figure  8.  Reflection  patterns  for  a  source  in  the  Kunlun,  northern  Tibet  at  phase  velocities 
of  3.7  and  4.3  km/s.  Rays  south  of  the  Himalayas  in  (a)  are  characterized  by  increased 
phase  velocities  and  widespread  transmission  into  the  mantle  due  to  constriction  of  the 
crustal  wave  guide.  The  reflection  pattern  in  (b)  is  highly  distorted;  both  transmission  into 
the  mantle  and  reflection  along  the  margins  of  Tibet  contribute  to  very  little  coherent  energy 
actually  exiting  through  the  crustal  waveguide. 
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FIGURE  8  (b) 
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ABSTRACT 

The  tracking  of  guided  wave  trains  by  ray  tracing  in  3-D  structures  can  provide  a  means 
of  interpreting  complex  seismograms  at  regional  distances.  The  method  relies  on  the 
interpretation  of  Lg  as  the  constructive  interference  of  multiple  S  reflections  within  the 
crust,  and  strong  scattering  can  be  simulated  by  the  inclusion  of  secondary  sources. 

This  approach  has  been  applied  to  Californian  events  observed  in  the  southwestern 
United  States  whose  records  at  regional  distances  frequently  exhibit  an  extended  and 
complex  Lg  coda  with  significant  late  energy  arriving  with  group  velocities  of  2.5  km/s 
or  less.  The  general  character  of  this  energy  precludes  an  interpretation  as  the  result  of 
stochastic  scattering  processes  or  dispersion  in  low  velocity  surface  sediments;  and  so 
some  alternative  explanation  needs  to  be  sought. 

The  observed  wavetrains  of  events  from  two  sites  in  southwestern  California 
(Coalinga  and  North  Baja)  at  stations  LAC,  MNV,  ELK  have  been  compared  with  the 
predictions  from  tracking  guided  wave  patterns.  Multipathing  of  energy  by  major 
changes  in  the  thickness  of  the  crustal  wave  guide  is  consistent  with  the  extended  and 
often  "pulse-like"  nature  of  the  Lg  coda.  Examples  of  plausible  multipathing 
mechanisms  include  reflection  from  the  ocean-continent  transition,  scattering  associated 
with  topographic  features  such  as  the  Sierra  Nevada,  and  resonance  within  the  narrow- 
corridor  of  Baja  California. 
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INTRODUCTION 

Bostock  &  Kennett  (1990)  have  recently  introduced  a  simple  scheme  for  modelling  the 
main  characteristics  of  guided  wave  trains  in  three-dimensionally  varying  structures  by 
tracking  the  behaviour  of  multiply  reflected  waves.  This  technique  has  been  able  to 
match  the  general  character  of  Lg  wave  propagation  in  central  Asia.  We  show  here  how 
the  method  can  be  adapted  to  examine  the  nature  of  the  propagation  paths  which  may 
contribute  to  the  L.g  wave  coda  in  the  southwestern  United  States,  by  including 
secondary  sources  to  simulate  scattering  as  well  as  guided  wave  propagation. 

The  stations  of  the  digital  seismic  network  operated  by  the  Lawrence  Livermore 
National  Laboratory  lie  at  regional  ranges  from  earthquakes  in  California  and  the 
northern  part  of  Baja  California.  The  most  prominent  part  of  the  seismic  record  is 
usually  the  Lg  phase  (with  frequencies  of  1  Hz  or  above),  although  lower  frequency  Rg 
waves  (around  0.25  Hz)  are  sometimes  the  largest  feature  (Xie  &  Mitchell  1990).  For  a 
number  of  events  the  structure  of  the  Lg  coda  is  quite  complex  with  distinct  packets  of 
energy  which  arrive  at  rather  low  apparent  group  velocities  (less  than  2.5  km/s).  In 
deriving  coda  Q  models  for  the  southwestern  United  States,  Xie  &c  Mitchell  (1990) 
excluded  such  records  from  their  analysis  since  this  class  of  arrivals  cannot  be  explained 
by  a  stochastic  scattering  model.  We  therefore  seek  a  more  direct  physical  interpretation. 

We  have  built  an  approximate  model  of  the  crustal  structure  in  California  and  Nevada 
based  on  the  5'  by  5'  digital  topography  database  (ET0P05).  The  shape  of  the  Moho  was 
derived  from  the  surface  topography  by  assuming  that  there  was  locally  approximately  80 
per  cent  isostatic  compensation  and  that  the  Moho  surface  was  smooth.  This 
compensation  factor  reproduces  the  Moho  depths  of  30-35  km  found  by  Hadley  & 
Kanamori  (1977)  beneath  the  Transverse  Ranges  of  southern  California,  but  may 
overestimate  crustal  thicknesses  for  offshore  structures.  Since  we  wish  to  reproduce  the 
main  features  of  the  propagation  we  do  not  attempt  to  make  a  full  allowance  for  the 
seismic  velocity  distribution  but  rather  impose  a  single  constant  velocity  of  3.5  km/s  in 
the  crust  and  a  velocity  of  4.6  km/s  beneath  the  crust-mantle  boundary.  We  then 
concentrate  on  the  way  in  which  the  guided  S  waves  are  affected  by  the  gradients  in 
surface  and  Moho  topography. 

The  character  of  the  guided  wavefield  for  different  phase  velocities  can  be  predicted  by 
changing  the  take-off  angles  of  the  rays  at  the  source.  The  propagation  pattern  is  then 
most  easily  followed  by  mapping  the  reflection  points  at  the  surface  and  at  the  base  of  the 
crust.  We  also  calculate  the  transit  times  along  each  ray  so  that  group  velocities  can  be 
calculated  along  the  propagation  paths.  We  consider  here  two  groups  of  earthquakes. 
7he  first  group  of  events  is  a  cluster  near  Coalinga,  in  the  Diablo  Range  of  central 
California,  recorded  at  stations  MNV  (315  km)  and  F,LK  (675  km)  in  Nevada.  The  records 
at  these  stations  show  a  complex  elongated  coda  with  a  few  distinct  later  arrivals  wliich 
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may  be  associated  with  back  reflection  from  the  continent-ocean  transition.  The  second 
group  of  events  lies  in  northern  Baja  California  recorded  at  station  LAC  (300  km)  and 
MNV  (775  km).  The  relatively  complex  records  suggest  significant  multipathing  and 
possibly  some  sort  of  resonator  effect  within  the  peninsula  of  Baja  California.  The 
locations  of  the  sources  and  receivers  are  given  in  Table  1  and  illustrated  in  fig.  1 
superimposed  on  the  topography  of  the  southwestern  United  States. 

DIABLO  RANGE  EVENTS: 

The  character  of  the  events  from  near  Coalinga  in  the  Diablo  Range  of  central 
California  is  illustrated  by  the  records  at  MNV  and  ELK  of  three  events  in  1983  with 
similar  epicentres  (Table  1).  These  displays  show  300  seconds  of  bandpassed 
seismograms  (corners  0.5,  1.5  Hz)  which  cover  the  Lg  phase  and  coda,  this  frequency 
window  is  designed  to  emphasize  Lg  type  energy.  The  three  events  have  significant 
variation  in  frequency  content  which  is  reflected  in  the  different  waveforms  for  Lg  in  the 
band  passed  seismograms.  The  times  corresponding  to  group  velocities  of  3.5  and  2.5 
km/s  are  indicated  by  markers.  We  note  that  there  is  still  substantial  energy  arriving 
very  late  in  the  coda. 

The  two  stations  lie  along  the  same  great  circle  and  the  propagation  path  crosses  the 
San  Joaquin  Valley  and  the  Sierra  Nevada.  There  is  the  possibility  of  some  energy  being 
delayed  by  trapping  in  the  sediments  of  the  San  Joaquin  Valley  but  this  would  not  be 
expected  to  produce  late  discrete  pulses  such  as  are  seen  at  both  MNV  and  ELK.  It  is  more 
likely  that  such  pulses  would  be  produced  by  scattering  or  reflection  mechanisms  giving 
rise  to  arrivals  reaching  the  stations  along  many  different  paths. 

The  great  circle  paths  for  the  three  Coalinga  events  are  very  close,  and  so  the  records 
should  reflect  similar  propagation  processes.  However,  the  significance  of  different 
processes  can  be  modified  by  the  character  of  the  radiation  pattern  for  the  event,  which 
will  be  affected  by  both  the  source  mechanism  and  depth. 

We  have  investigated  the  nature  of  the  propagation  processes  which  might  contribute 
to  the  character  of  the  observed  records  using  the  ray-based  method  of  Bostock  &  Kennetl 
(1990).  In  figure  3  we  show  the  propagation  pattern,  from  the  source,  to  the  two  receivers 
for  both  direct  arriving  waves  and  waves  reflected  back  from  the  continent-ocean 
transition  just  off-shore.  The  phase  velocity  has  been  chosen  to  be  4.0  km/s  to  give  a 
good  general  representation  of  the  propagation  behaviour.  The  ray  paths  through  the 
model  are  shown  and  superimposed  on  these  are  small  squares  marking  the  reflection 
points  for  the  multiply  reflected  S  waves  at  the  base  of  the  crust.  The  squares  are  filled 
where  reflection  is  subcritical  to  indicate  the  possibility  of  energy  escaping  from  the 
waveguide. 
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In  the  case  of  direct  propagation,  the  wavefronts  marked  by  the  patterns  of  reflection 
points  are  coherent  for  the  first  two  crustal  reflections  but  become  more  distorted  after 
surface  reflections  have  occurred  within  the  Sierra  Nevada.  At  station  ELK  the  pattern 
of  reflections  for  nearby  rays  shows  considerable  variation  and  this  would  lead  to  a 
spread  in  the  apparent  group  velocities  for  the  main  energy. 

The  Coalinga  events  are  not  very  far  from  the  coast  and  the  sharp  contrast  between 
continental  and  oceanic  structures  can  lead  to  reflection  of  the  guided  waves.  Because  of 
the  way  we  have  built  our  crustal  model  the  details  of  the  transition  may  be  somewhat 
simplified  but  the  general  character  of  the  reflection  process  can  be  quite  well 
represented.  Energy  can  be  returned  from  the  coastal  zone  to  both  MNV  and  ELK  from  a 
number  of  different  locations.  We  have  considered  the  full  range  of  phase  velocities 
which  might  be  regarded  as  of  significance  for  the  Lg  waves  and  we  find  that  the  reflected 
arrivals  would  be  expected  to  arrive  at  MNV  with  apparent  group  velocities  from  1.8 
km/s  to  1.6  km/s.  For  ELK  the  corresponding  group  velocity  range  is  from  2.3  km/s  to 
2.2  km/s  (since  the  station  is  further  away  the  additional  path  length  to  the  reflection 
point  and  back  is  a  smaller  fraction  of  the  whole).  The  spread  of  times  predicted  for  such 
reflections  is  indicated  on  each  of  the  seismograms  in  fig.  2  in  a  grey  bar.  For  both  MNV 
and  ELK  there  is  an  interruption  of  the  decay  pattern  for  the  Lg  coda  in  the  predicted 
time  span  on  the  band-passed  records,  with,  in  most  cases,  a  distinct  increase  in 
amplitude.  It  is  difficult  to  model  the  strength  of  such  a  reflection  from  the 
continent/ocean  transition  in  detail  because  of  strong  dependence  on  the  character  of  the 
source  radiation,  but  the  reflectivity  is  unlikely  to  exceed  0.20  at  the  angles  of  incidence 
in  fig.  3.  The  reflection  would  therefore  only  be  expected  to  be  recognized  when  the 
general  amplitude  level  is  low.  At  MNV  the  predicted  time  interval  on  the  unfiltered 
records  is  masked  by  low  frequency  Rg.  However  at  ELK  the  interval  lies  in  the  coda  of 
Rg  and  there  is  a  distinct  arrival  at  this  time  (fig.  4)  on  the  unfiltered  records,  w'hich  is 
most  pronounced  for  the  event  of  1983  July  9. 

We  have  already  noted  that  the  pattern  of  the  guided  waves  becomes  disrupted  on 
passage  through  the  Sierra  Nevada.  The  complex  topography  and  structure  in  this  zone 
will  also  generate  scattered  waves  which  would  tend  to  extend  the  coda  of  the  guided 
waves.  We  have  introduced  a  simple  model  of  a  number  of  fixed  secondary  scatterers 
along  the  Sierras  and  then  look  at  the  time  relations  and  propagation  paths  required  to 
get  energy  from  the  source  to  the  receivers  via  the  secondary  scattering  centres.  The 
process  is  illustrated  in  figure  4  for  the  Coalinga  source  and  a  receiver  at  ELK.  For  waves 
travelling  with  a  true  group  velocity  of  3.5  km/s,  the  scattering  process  can  easily  extend 
the  train  to  3.0  km/s.  For  MNV  the  contribution  of  the  scattered  wave  legs  is  much 
larger  and  the  apparent  group  velocity  can  be  reduced  to  2.5  km/s.  Similar  stretch  factors 
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will  occur  for  other  group  velocities.  The  corrugated  surface  topography  in  the  Basin 
and  Range  province  is  also  likely  to  contribute  to  scattering  processes. 

As  has  often  been  suggested  (see  e.g.  Aki  (982),  Herraiz  &  Espinosa  (1987))  scattering  is  a 
powerful  mechanism  for  fueling  energy  into  the  latter  parts  of  the  seismogram. 
Stochastic  modelling  of  such  scattering  depends  on  an  ensemble  average  over  many 
models  and  can  therefore  predict  envelope  decay  rates  averaged  over  many  seismograms 
but  not  the  details  of  a  single  record  (Hudson  1982).  More  specific  modelling  is  therefore 
needed  to  try  to  unravel  the  nature  of  the  Lg  coda. 

NORTHERN  BAJA  CALIFORNIA  EVENTS 

The  second  group  of  events  occurred  in  northern  Baja  California  and  were  recorded  at 
the  stations  LAC  (near  Palm  Springs)  and  MNV  which  lie  close  to  a  great  circle  through 
the  event.  Whereas  the  Coalinga  events  showed  fairly  consistent  records,  there  is  more 
variability  from  this  group  of  events.  At  LAC,  in  particular,  the  event  character  seems  to 
vary  between  relative  rapid  decay  and  a  more  extended  coda.  But  when  examined  in 
detail  the  similarity  between  different  events  of  the  second  group  is  somewhat  greater; 
the  coda  includes  short  bursts  of  energy  which  continue  with  diminished  amplitude  for 
150  seconds. 

In  figure  6  we  show  300  s  of  bandpassed  seismograms  for  the  Lg  wave  and  its  coda  for 
the  event  of  1985  May  10  at  MNV  and  LAC,  with  again  the  times  for  apparent  group 
velocities  of  3.5  km/s  and  2.5  km/s  marked.  The  striking  features  of  the  records  at  MNV 
are  the  complex  nature  of  the  main  arrivals  and  the  persistence  of  the  coda  to  very  long 
times.  Similarly  at  LAC,  although  the  main  coda  decay  is  quite  rapid  the  level  remains 
sustained  at  a  moderate  level  for  nearly  50  s  (from  200  to  250  s)  before  declining  further. 

The  application  of  the  Bostock  &  Kennett  technique  to  this  class  of  event  showed  that 
the  propagation  patterns  were  quite  sensitive  to  the  location  of  the  event,  which  would 
fit  with  the  variations  in  the  character  of  the  seismograms  from  different  events.  We 
illustrate  in  figure  7  the  propagation  pattern  for  the  direct  arrivals  and  for  waves 
reflected  at  grazing  incidence  from  the  continent/ocean  transition,  as  well  as  a  few  waves 
that  get  trapped  in  the  narrow  corridor  of  Baja  California. 

The  spine  of  the  Sierras  tends  to  trap  waves  travelling  in  the  general  direction  of 
MNV,  these  waves  can  then  scatter  to  give  many  possible  paths  by  which  energy  can 
arrive  at  MNV  within  the  group  velocity  window  from  3.5  km/s  to  3.0  km/s.  In 
addition  there  is  a  group  of  ray  paths  to  MNV  corresponding  to  reflection  from  offshore 
in  the  neighbourhood  of  the  Channel  Islands;  such  waves  would  arrive  w'ith  an 
apparent  group  velocity  around  2. 4-2. 3  km/s.  3  he  predicted  time  interval  for  these 
reflections  is  marked  in  grey  on  figure  5  and  does  indeed  correspond  to  a  distinct  later 
arrival  at  MNV.  For  a  single  event  such  agreement  might  be  accidental.  However  in  fig. 
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8  unfiltered  records  at  MNV  for  two  Baja  events  with  similar  locations  are  shown,  both 
have  distinct  arrivals  in  the  predicted  lime  interval.  Since  the  reflection  is  near-grazing, 
quite  high  reflection  amplitudes  would  be  feasible.  For  station  LAC  the  required  angle  of 
reflection  from  the  offshore  structure  is  too  acute  to  enable  a  consistent  pattern  of  rays  to 
be  established,  and  so  pulse  like  arrivals  would  not  be  expected. 

The  sources  lie  at  the  northern  end  of  the  narrow  peninsula  of  Baja  California  and  it  is 
tempting  to  suggest  that  the  persistence  of  coda  energy  on  the  records  is  due  to  some 
kind  of  resonance  phenomenon  akin  to  a  organ-pipe  mode  within  the  peninsula 
returning  energy  over  a  long  time  span  as  can  be  clearly  seen  in  fig.  8.  Ray  tracing 
calculations  (as  in  figure  7)  confirm  that  it  is  indeed  possible  to  trap  energy  within  Baja 
California  but  that  the  configuration  depends  strongly  on  the  assumptions  used  in 
building  the  crustal  model.  In  addition  to  the  rays  presented  in  this  study  some 
scattering  may  occur  where  the  trend  of  the  structure  alters  and  such  waves  would  arrive 
at  about  the  right  time  to  give  the  extended  main  arrivals  seen  in  figure  5. 

DISCUSSION 

We  have  shown  that  it  is  possible  to  explain  much  of  the  character  of  observed 
regional  S  wave  seismograms  in  a  complex  region  such  as  the  southwestern  United 
States  using  calculations  based  on  Bostock  &  Kennett's  (1990)  treatment  of  guided-wave 
propagation  in  three-dimensionally  varying  media  for  a  simple  crustal  model.  The 
inclusion  of  scattering  and  reflection  phenomena  as  well  as  direct  propagation  help  to 
provide  the  classes  of  paths  needed  to  span  the  observed  ranges  of  group  velocity.  The 
continent/ocean  transition  is  a  very  effective  barrier  to  Lg  wave  propagation  and  can 
give  rise  to  quite  efficient  reflections.  Mountain  ranges  tend  to  encourage  ducted  wave 
propagation  and  where  such  a  feature  is  coupled  with  efficient  reflecting  margins  such  as 
in  Baja  California  it  is  possible  to  build  up  a  strong  trapped  wavefield,  w'hich  may  act  as  a 
resonator  feeding  energy  back  into  the  external  seismic  waves  over  a  fairly  long  time 
span. 

The  guided  wave  tracking  procedure  can  be  used  to  predict  the  character  of  regional 
phases  given  a  knowledge  of  crustal  structure.  Indeed  such  modelling  had  suggested  the 
resonator  character  of  Baja  California  before  any  observations  were  sought  to  check  the 
result. 
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TABLE  1.  Locations  of  seismic  stations  and  events  studied 


Location  ofLLNL  Stations: 


ELK 

115.2388°  W 

40.7448°  N 

MNV 

118.1544°  W 

38.4322°  N 

LAC 

116.4115°  W 

34.3898°  N 

Events: 

Group  1: 

Coalinga,  Diablo  Range 

83  vii  09 

7:40:50.90 

120.41°  W 

36.24°  N 

depth  10  km 

mb  5.3 

83  vii  22 

3:43:00.60 

120.41°  W 

36.21°  N 

depth  10  km 

mb  5.3 

83  ix  1 1 

11:48:08.00 

120.39°  W 

36.23°  N 

depth 

8  km 

mb  5.0 

Group  2: 

Northern  Baja  California 

85  V  08 

23:40:18.20 

115.80°  W 

31.70°  N 

depth 

5  km 

mb  5.1 

85  V  10 

5:39:36.30 

115.77°  W 

31.70°  N 

depth 

6  km 

mb  4.6 
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FIGURE  CAPTIONS; 


1.  Reference  map  of  southwestern  United  States  with  locations  of  the  stations  of  the  LLNL  network 
and  events  uscxi  in  this  study. 

2.  Bandpassed  seismograms  of  the  Lg  wave  and  coda  for  the  group  of  events  in  the  coastal  range  of 
central  California  at  the  stations  MNV  and  ELK  in  Nevada.  The  passband  (0.5,  1.5  Hz)  is 
designed  to  enhance  the  Lg  wave.  The  time  corresponding  to  a  group  velocity  of  3.5  km/s  is 
marked  by  a  square  and  2.5  km/s  by  a  triangle.  The  predicted  time  span  for  reflections  from  the 
continent  ocean  transition  are  marked  by  a  grey  bar. 

3.  Propagation  diagram  for  guided  S  waves  from  a  source  near  Coalinga  at  a  phase  velocity  of  4.0 
km/s.  Both  direct  propagation  and  reflection  from  the  continent/ocean  transition  are  illustrated. 
The  location  of  MNV  and  ELK  are  shown  by  solid  triangles  and  the  source  by  a  hexagon, 

4  Unfiltered  seismograms  for  the  group  of  Coalinga  events  at  ELK.  The  predicted  tiinc  span  for 
reflections  from  the  continent/ocean  transition  arc  marked  by  a  grey  bar  and  correspond  to 
distinct  arrivals  in  the  coda. 

5.  Simulation  of  scattering  by  inclusion  of  secondary  sources  along  the  Sierra  Nevada,  the 
propagation  patterns  shown  arc  for  a  receiver  at  ELK. 

6.  Bandpassed  seismograms  of  the  Lg  wave  and  coda  for  the  event  of  1985  May  8  in  the  northern  Baja 
California  at  the  stations  LAC  (near  Palm  Springs)  and  MNV  in  Nevada.  The  passband  (0.5,  1.5 
Hz)  is  designed  to  enhance  the  Lg  wave.  The  time  corresponding  to  a  group  velocity  of  3.5  km/s  is 
marked  by  a  square  and  2.5  km/s  by  a  triangle.  The  predicted  time  span  for  reflections  from  the 
continent  ocean  transition  are  marked  by  a  grey  bar. 

7.  Propagation  diagram  for  guided  S  waves  from  a  source  in  northern  Baja  at  a  phase  velocity  of  4.0 
km/s.  Direct  propagation  and  grazing  reflection  from  the  continent/ocean  transition  are 
illustrated,  as  well  as  the  propensity  for  rays  to  get  trapped  within  Baja  California.  The 
location  of  LAC  and  MNV  are  shown  by  solid  triangles  and  the  source  by  a  hexagon. 

8.  Unfiltered  seismograms  for  two  northern  Baja  events  at  MNV.  The  predicted  time  span  for 
reflections  from  the  continent/ocean  transition  in  the  neighbourhood  of  the  Channel  Islands  are 
marked  by  a  grey  bar  and  correspond  to  distinct  arrivals  in  the  coda. 
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THE  EFFECTS  OF  CRUSTAL  VELOCITY  GRADIENTS  ON  THE  PROPAGATION  OF 
REGIONAL  S  PHASES 


].  Roger  Bowman  and  B.L.N.  Kennett 
Research  School  of  Earth  Sciences, 
Australian  National  University 
GPO  Box  4,  Canberra,  ACT  2601 
Australia 


INTRODUCTION 

Lg  waves  are  the  dominant  feature  of  seismograms  at  regional  distances  and  are 
therefore  quite  useful  for  estimating  earthquake  magnitudes  and  for  the  detection 
and  estimation  of  the  yields  of  nuclear  explosions.  However,  the  accuracy  of 
magnitude  and  yield  estimates  is  dependent  on  the  correct  estimation  of  Lg 
amplitude  decay  with  distance.  In  this  paper  we  present  results  of  a  study  of  Lg 
attenuation  in  the  Precambrian  shield  in  Central  Australia  that  demonstrate  the 
effect  of  regional  crustal  structure  on  the  rate  of  Lg  amplitude  decay. 

Aftershocks  of  the  22  January  1988  Tennant  Creek,  Northern  Territory 
earthquakes  provide  a  localized  source  of  radiation  for  an  analysis  of  regional 
phases  recorded  by  a  425  km  linear  array  of  portable  seismographs  between 
Tennant  Creek  and  Alice  Springs  (Figure  1).  The  three  Tennant  Creek  mainshocks 
(Ms  6.3,  6.4  and  6.7)  produced  two  scarps  of  32  km  combined  length  (Bowman,  1988) 
and  were  followed  by  thousands  of  aftershocks  in  a  10  km  x  40  km  zone  (Bowman 
et  al.,  1990)  and  in  a  complex  geometry  that  suggests  faulting  on  conjugate  planes 
(Choy  and  Bowman,  1990). 

Most  previous  studies  of  Lg  propagation  have  relied  on  broadly  distributed 
sources  and  receivers  and,  as  a  consequence,  estimates  of  phase  velocity  and 
attenuation  represent  averages  over  the  region  sampled  by  the  paths  available. 
Moreover,  there  is  seldom  redundant  data  for  the  same  paths,  so  no  direct  measure 
of  errors  in  the  estimates  is  available.  We  take  advantage  of  the  abundant 
aftershocks  of  the  Tennant  Creek  earthquakes  to  examine  Lg  propagation  along  a 
simple  path  in  the  North  Australian  Craton. 

Several  lines  of  evidence  suggest  that  the  North  Australian  craton  has  low 
seismic  attenuation.  There  is  a  large  positive  bias  in  Ml  estimates  for  regional  and 
teleseismic  earthquakes  recorded  at  stations  located  on  the  Australian  craton  that 
may  be  attributed  to  low  attenuation  (e  g.,  Gaull  et  al.,  1989).  Also,  deep  seismic 
reflection  data  from  central  Australia  is  unusually  rich  in  high  frequencies  (Goleby 
et  al.,  1989).  Furthermore,  the  isoseismal  areas  of  earthquakes  in  the  Australian 
craton  are  comparable  to  those  in  other  stable  continental  regions  with  low 
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attenuation.  However,  in  the  present  study,  we  find  that  estimates  of  attenuation 
are  anomalously  high  when  we  correct  the  observed  Lg  wave  amplitudes  for 
geometrical  spreading  with  a  factor  commonly  used  in  North  America  and 
Western  Europe.  Crustal  studies  in  central  Australia  suggest  that  the  crust-mantle 
transition  consists  of  a  gradient  zone  rather  than  an  abrupt  discontinuity 
(Finlayson,  1982).  We  have  therefore  investigated  the  effect  of  the  character  of  the 
crust-mantle  transition  on  the  apparent  attenuation  of  Lg  waves  using  synthetic 
seismogram  modelling.  Incomplete  trapping  of  Lg  waves  in  the  presence  of  a 
gradient  zone  is  likely  to  be  responsible  for  the  relatively  rapid  decrease  in  Lg 
amplitude  with  distance  from  the  source. 

DATA  AND  GEOLCXIICAL  SETTING 

Eight  portable  seismographs  were  operated  on  a  profile  between  Tennant  Creek 
and  Alice  Springs  for  one  week  in  March,  1988  (Figure  1).  These  stations  consisted 
of  vertical-component  Willmore  Mark  111  seismometers  with  1  s  natural  periods 
and  9-track  digital  instruments  recording  continuously  at  25.6  samples/s.  A 
tripartite  array  consisting  of  station  N1  situated  west  of  the  fault  zone,  three- 
component  station  B2  (operated  at  16  samples/s)  at  the  Warramunga  (WRA)  array 
and  the  northernmost  station  HI  of  the  linear  array  was  used  for  determination  of 
aftershock  epicenters.  Data  from  station  H2,  which  had  a  partially  clamped 
seismometer,  and  station  H4,  which  had  numerous  tape  errors,  were  not  used  for 
the  attenuation  analysis,  but  data  from  H2  were  retained  in  the  record  section 
displays.  In  addition  to  the  portable  stations,  we  use  data  from  the  three- 
component,  short-period  station  (SZ,  SN  and  SE  for  vertical,  north-south  and  east- 
west  components)  of  the  Alice  Springs  array  (A5PA,  Figure  1). 

The  propagation  path  from  Tennant  Creek  to  Alice  Springs  is  in  the  North 
Australian  Craton,  which  has  formed  a  single  relatively  stable  crustal  block  for  1700 
m.y.  (Plumb,  1979).  The  Tennant  Creek  fault  zone  is  near  the  boundary  between 
the  Paleozoic  Wiso  Basin  and  the  Early  Proterozoic  Tennant  Creek  "inlier".  The 
inlier  is  a  basement  outcrop  of  a  much  larger  feature  and  consists  of  sediments  and 
metasediments  of  the  Warramunga  Group,  which  are  intruded  by  lower 
Proterozoic  granites.  Stations  H1-H3  were  located  on  the  Proterozoic  Davenport 
geosyncline,  stations  H5  and  H6  were  situated  on  the  western  flank  of  the  late 
Proterozoic-Devonian  Georgina  Basin  and  stations  H7,  H8  and  ASPA  were  on  the 
Proterozoic  Arunta  block.  None  of  the  stations  is  underlain  by  more  than  50-100  m 
of  sediment.  All  seismometers  were  emplaced  in  the  Cenozoic  aeolean  and 
residual  sand  and  gravel  deposits  that  blanket  most  of  the  region  with  a  thickness 
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generally  less  than  20  in.  Therefore,  there  is  no  reason  to  expect  large  variations  in 
site  response. 

Because  of  the  high  level  of  aftershock  activity,  we  were  able  to  select  11 
earthquakes  for  analysis  of  sufficient  size  to  be  recorded  at  the  most  distant  stations 
during  the  brief  recording  period  (Figure  2,  Table  1).  Aftershock  locations  were 
determined  using  arrival  time  data  from  stations  Nl,  B2  and  HI,  a  simplified 
model  (Table  2)  based  on  refraction  surveys  (Finlayson,  1981,  1982),  and  the 
program  HYPOELLIPSE  (Lahr,  1980).  Additional  arrival  times  were  available  from 
triggered  recorders  in  the  source  area  for  several  of  the  larger  events.  Because  only 
three  stations  were  generally  used  in  the  solution  and  no  station  was  located 
within  two  or  three  focal  depths,  all  but  one  of  the  depths  were  fixed  at  4  km, 
which  is  half  the  maximum  focal  depth  observed  in  detailed  aftershock  studies 
(Bowman  et  al.  ,  1990).  Uncertainties  in  epicenter  are  generally  less  than  3  km. 

For  event  5  (Table  1)  a  sample  record  section  of  the  regional  data  is  shown  in 
Figure  3  and  a  three-component  record  from  Alice  Springs  is  shown  in  Figure  4.  In 
both  displays  each  trace  is  normalized  to  its  maximum  amplitude.  At  Alice 
Springs  the  north-south  component  (SN)  is  very  close  to  radial  and  the  east-west 
component  (SE),  to  transverse.  The  vertical  component  data  in  Figure  3  shows  the 
evolution  of  the  wavefield  from  10  km  to  430  km  from  the  source.  At  ranges 
beyond  120  km,  the  seismograms  consist  of  impulsive  crustal  P  waves  followed  by 
emergent  Lg  wavetrains.  The  Lg  arrival  at  station  H6  is  more  impulsive  and 
stronger  than  at  neighboring  stations.  Record  sections  for  the  11  events  selected  are 
remarkably  similar,  although  there  is  some  variation  in  the  relative  strengths  of  P 
and  Lg,  probably  caused  by  differences  in  focal  mechanism.  The  smaller 
earthquakes  have  the  highest  signal-to-noise  ratio  in  the  band  1.5-5  Hz. 

The  P-wave  first  arrivals  have  an  apparent  velocity  of  7. 7-7. 8  km/s,  which 
suggests  that  they  are  refractions  from  the  lower  crust.  Finlayson  (1982)  used  data 
from  mine  blasts  near  Tennant  Creek  recorded  along  a  line  extending  500  km  to 
the  east  to  interpret  a  Moho  depth  of  50  km  and  Pn  velocity  of  8.2  km/s.  The 
crossover  distance  where  Pn  becomes  a  first  arrival  in  Finlayson’s  model  TCMI-2  is 
about  350  km  for  a  surface  source,  and  only  stations  H8  and  ASPA  lie  beyond  this 
distance.  However,  the  differences  between  arrival  times  at  these  stations  in  our 
data  are  consistent  with  the  velocities  of  7.7-7.8  km/s  observed  at  closer  ranges. 
This  suggests  that  the  crust  is  thicker  along  this  profile  extending  south  of  Tennant 
Creek  than  along  a  profile  to  the  east,  or  that  the  Pn  velocity  is  anomalously  low. 
The  apparent  velocity  of  7.7-7. 8  km/s  observed  here  is  intermediate  betw'een  the 
lower  crustal  (7.4  km/s)  and  Pn  (8.2  km/s)  velocities  of  TCMI-2. 
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On  the  vertical  components  along  the  profile  (Figure  3),  the  initial  P  wave  is 
followed  by  a  coda  of  rather  uniform  amplitude  without  prominent,  discrete 
phases  such  as  Pg  or  Sn  (except  for  a  clear  early  Sn  arrival  at  H8).  On  the  radial 
(SN)  component  at  ASPA  (Figure  4),  on  the  other  hand,  the  first  P  is  followed  by  »1 
s  by  a  Pg  phase,  and  on  both  the  horizontal  components,  the  Lg  is  preceded  by  a 
smaller  pulse  that  may  be  an  S  counterpart  to  the  initial  P  onset  arrival.  It  should 
be  noted  that  trace  normalization  exaggerates  the  strength  of  the  Pg  on  the  radial 
(SN)  component. 

The  Lg  phase  dominates  the  seismograms  at  distances  beyond  120  km  with  an 
average  amplitude  3.3±1.4  times  that  of  the  P  wave  in  a  group  velocity  window 
from  8-7  km/s  (for  seismograms  filtered  in  a  3-8  Hz  passband).  Unless  otherwise 
specified,  amplitudes  refer  to  the  maximum  sustained  amplitude  defined  by  Nuttli 
(1973),  that  is  the  amplitude  equaled  or  exceeded  by  the  three  largest  amplitude 
peaks  in  the  wavetrain.  The  Lg/P  ratio  in  northern  Australia  falls  between  that 
observed  for  earthquakes  in  eastern  North  America  (10)  and  for  earthquakes  and 
explosions  in  the  western  USSR  (1)  (Pomeroy  et  al.,  1982).  At  ASPA  the  ratios  of 
radial  and  transverse  components  to  the  vertical  component  amplitudes  are 
0.8±0.2  and  1.3±0.3,  respectively,  somewhat  lower  than  the  ratio  of  2  between 
maximum  horizontal  and  vertical  sometimes  assumed  (e  g..  Street,  1978). 

GROUP  VELOCITY 

The  Lg  group  velocity  is  estimated  as  3.7  km/s  by  fitting  a  line  to  the  maximum 
amplitudes  recorded  at  stations  H3  to  ASPA  for  the  1 1  earthquakes  in  Table  1,  and  a 
similar  velocity  is  estimated  from  the  onsets  of  the  Lg  wavetrains.  The  group 
velocity  in  Central  Australia  is  at  the  upper  end  of  the  range  usually  reported  for 
stable  continents  (Pomeroy  et  al.,  1982),  and  is  somewhat  higher  than  the  only 
previous  estimate  of  3.5  km/s  for  Australia,  based  on  paths  crossing  both 
Precambrian  shield  and  the  Phanerozoic  region  of  eastern  Australia  (Bolt,  1957). 

LG  ATTENUATION 

The  Lg  amplitude  at  frequency  /  for  the  ilh  earthquake  (/=l,n)  observed  at  the  /th 
station  (;=l,w)  can  be  approximated  by 

A,,(f)  =  Ao,(f)  i;(/)  R,<f)  S,,  D,,  e- r,,,  ( 1 ) 

where  Ad  is  a  source  factor,  Ij  is  the  instrument  response,  Rj  is  the  site  response,  Djj 
is  a  geometric  spreading  factor,  Sij  is  a  source  radiation  term  and  e'/’'i/  is  the 
attenuation  factor  for  source-station  separation  r,j.  The  coefficient  of  attenuation 
Y=  nf  /QU,  where  Q  is  the  temporal  quality  factor  and  U  the  group  velocity,  and  rij 
is  the  source-receiver  distance.  If  the  instruments  are  nominally  identical 
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log  Aij  -  log  Dij  -  log  Sij  =  log  Aoi  -dogioe)  yr,;.  (2) 

where  the  site  response  has  been  taken  as  unity.  Now,  define  Anj*  as  the 
logarithmic  amplitude  corrected  for  geometrical  spreading  and  source  radiation 
effects,  and  set  y'  =  (logige)  y.  Under  the  assumption  that  a  single  attenuation 
coefficient  is  adequate  to  describe  the  nearly  coincident  paths  on  the  Tennant  Creek 
profile,  and  variations  in  source  radiation  patterns  can  be  ignored. 

Ai*  =  \ogAoi-r'rij  (3) 

With  a  number  of  amplitude  observations  we  can  estimate  the  spatial  attenuation 
coefficient  y'  and  the  source  excitation  factors  Ad  by  a  generalised  least  squares 

technique.  The  estimates  are  determined  by  minimising 

m  n 

=  I  I  Wij  (/lJ-Mo.-y'r.;))2  (4) 

where  Wij  is  a  weighting  factor  defined  to  reduce  the  influence  of  small  arrivals 
(Dwyer  et  al  1983) 

r  1  SNR>4 

Wij  =  ]  SNR/2-1  2<SNR<4  (5) 

[  0  SNR<2 

and  SNR  is  the  ratio  of  Lg  amplitude  to  the  noise  level  before  the  first  P  arrival. 

The  geometrical  spreading  factor  for  Lg  is  frequently  approximately  by  the  simple 
power  law 

D  =  r-n  (6) 

Nutti  (1973)  suggested  that  n  =  5/6.  This  is  the  predicted  amplitude  decay  with 

distance  for  an  Airy  phase  trapped  in  the  crustal  w-aveguide  (p  145,  Ewing  et  al 
1957).  However,  Lg  is  not  a  single  Airy  phase,  but  rather  a  superposition  of  Airy 
phases  of  different  modes  with  a  range  of  phase  and  group  velocities.  As  the 
observed  amplitude  will  depend  on  the  details  of  interference  of  the  constituent 
modes,  the  amplitude  decrease  may  not  follow  such  a  simple  relationship. 
Although,  results  of  numerical  experiments  by  Campillo  et  al.  (1984)  and  Shin  and 
Herrmann  (1987)  are  consistent  with  an  average  trend  close  to  r  ,  it  is 

noteworthy  that  both  of  these  studies  focused  on  velocity  models  with  sharp 
discontinuities  at  the  Moho.  Because  of  the  difficulty  in  simultaneously 
estimating  geometric  spreading  and  attenuation  owing  to  absorption  and 
scattering,  we  initially  assume  r  "5/6  spreading  and  arrive  at  an  unexpectedly  low- 
value  for  Q  in  Central  Australia. 

We  have  Lg  amplitudes  in  four  discrete  bands  with  center  frequencies  fc  of  1  Hz, 
2  Hz,  3  Hz  and  8  Hz  and  lower  and  upper  corner  frequencies  of  /c/'^2  and  V2  fc. 
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respectively.  For  the  smaller  events,  the  8  Hz  frequency  band  had  a  poor  SNR,  so 
we  concentrate  on  the  three  lower  frequency  bands.  The  solution  to  equation  (4) 
gives  an  intercept  for  each  earthquake,  proportional  to  the  size  of  the  event,  and  an 
estimate  of  y consistent  with  the  assumed  geometrical  spreading.  Maximum 
amplitude,  maximum  sustained  amplitude  and  RMS  amplitude  in  the  group- 
velocity  window  3.7-3. 1  km/s  were  used  as  measures  of  the  amplitude  and  found 
to  give  consistent  results.  An  example  of  the  least  squares  solution  for  the  1.41 
-2.82  Hz  passband  for  RMS  amplitudes  is  shown  in  Figure  5  and  implies  a  Q  of  302 
for  n=5/6. 

Such  a  relatively  low  Q  value  is  inconsistent  with  the  frequency  content  of  the 
seismograms,  and  so  it  is  appropriate  to  examine  the  validity  of  the  assumptions 
on  which  the  attenuation  estimate  is  based.  In  order  to  extract  the  loss  coefficient 
we  have  to  make  an  allowance  for  the  main  propagation  effects  of  wavefront 
spreading.  Figure  6  illustrates  the  trade-off  between  the  attenuation  coefficient  y 
and  the  index  n  in  the  power  law  representation  of  the  spreading  (6).  An 
appropriate  choice  of  spreading  relation  is  clearly  critical  for  an  accurate  estimate  of 
attenuation. 

In  the  following  section  we  show  how  the  pattern  of  Lg  amplitudes  variation 
with  distance  is  influenced  by  the  character  of  the  crust-mantle  transition.  A 
simple  power  law  relation  for  the  amplitude  variation  with  distance 
underestimates  the  complexity  of  the  Lg  propagation  process,  and  in  those  cases 
where  focussing  effects  give  localised  high  amplitudes  is  a  very  poor 
representation  of  the  behaviour. 

DISCUSSION 

As  noted  above,  crustal  models  for  the  Tennant  Creek  inlier  suggest  that  the 
crust-mantle  transition  is  a  gradient  zone  rather  than  an  abrupt  discontinuity 
(Finlayson  1981,  1982).  Estimates  of  the  spatial  attenuation  coefficient  for  Lg  are 
strongly  dependent  on  the  assumptions  made  in  correcting  for  the  major 
geometric  spreading  component.  The  choice  of  gives  a  noticeable  discrepancy 
between  the  low  estimated  Q  values  and  the  propagation  of  high  frequencies  to 
considerable  distances. 

In  order  to  investigate  the  influence  of  the  nature  of  the  crust-mantle  transition 
on  the  character  of  the  Lg  phase  we  have  constructed  synthetic  seismograms  for  a 
number  of  simple  models  of  the  crust-mantle  transition.  The  reference  model  SJU 
has  a  sharp  Moho  discontinuity  at  a  depth  of  38  km  and  two  classes  of  models  were 
constructed  based  on  replacing  the  single  discontinuity  by  a  set  of  three  steps. 
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In  the  sequence  SGI  to  SG3  (Figure  7a)  the  discontinuity  is  replaced  by  a  velocity 
gradient,  covering  the  full  span  of  the  original  Moho  jump,  extending  over  a 
progressively  larger  depth  interval. 

In  the  second  sequence  SJl,  SJ2,  SG3  (Figure  7b)  a  partial  jump  at  the  base  of  the 
crust  at  38  km  is  retained  and  the  remainder  of  the  velocity  increase  to  mantle 
velocities  is  in  the  form  of  a  gradient.  The  models  SJO  and  SG3  can  be  regarded  as 
the  limits  of  both  sequences. 

Theoretical  seismograms  were  calculated  using  a  wavenumber  integral  approach 
for  each  of  the  models  (Kennett  1983,  1988).  A  common  set  of  source  and 
integration  parameters  were  employed  to  allow  direct  comparison  between  the 
results.  The  calculations  were  carried  out  for  an  explosive  source  at  15  km  (  to 
avoid  source  radiation  effects),  with  a  P  wave  0  of  1000  and  an  S  wave  Q  of  500.  A 
frequency  interval  of  0.5-4  Hz  and  a  phase  velocity  interval  of  3.0-5.0  km/s  were 
used  in  the  integration,  and  all  P  and  S  wave  crustal  multiples  were  included.  The 
resulting  seismograms  for  the  two  sequences  of  models  are  displayed  in  Figure  8  in 
a  true  amplitude  display  with  a  reduction  velocity  of  4  km/s. 

At  the  head  of  each  panel  we  display  the  seismograms  for  the  reference  model 
SJO  with  the  sharp  Moho  transition.  Sn  is  seen  as  a  weak  arrival  at  200  km  with  a 
reduced  time  of  5  s  which  decreases  to  -2  s  at  400  km.  Sn  is  followed  by  a  similar 
arrival  with  an  extra  leg  reflected  above  the  source.  Such  a  free  surface  'ghost'  also 
helps  to  add  to  the  complexity  of  the  energy  propagating  in  the  crust.  The  sequence 
of  packets  of  multiply  reflected  S  wave  energy  within  the  crust  build  up  to  form  the 
Lg  wave  train.  The  interference  phenomena  leading  to  the  generation  of  Lg  do  not 
lead  to  a  monotonic  decay  of  amplitude  with  increasing  range.  The  largest 
amplitudes  occur  for  near  critical  reflection. 

In  the  sequence  of  gradient  models  (Figure  8a)  the  general  pattern  of  the 
seismograms  is  maintained  but  the  locus  of  the  maximum  amplitude  of  Lg  is 
shifted  to  greater  range  as  the  gradient  that  has  replaced  the  Moho  discontinuity  is 
reduced.  Because  of  the  reduced  gradient,  the  reflected  S  waves  from  the  crust- 
mantle  transition  spend  longer  in  this  region  and  so  require  a  greater  distance  to 
return  to  the  surface.  A  further  consequence  is  that  these  waves  travel  faster  than 
before  and  so  the  locus  of  maximum  amplitude  in  Lg  shifts  to  higher  group 
velocities. 

The  theoretical  seismograms  for  the  sequence  of  models  with  progressively 
smaller  velocity  jump  at  38  km  depth  are  displaced  in  Figure  8b.  In  model  SJl  the 
large  velocity  jump  at  38  km  means  that  critical  reflections  are  important  even 
though  Sn  is  enhanced  from  model  SJO  by  the  increased  gradient  below  38  km.  For 
model  SI?  with  a  smaller  velocity  jump  but  increased  gradient  below  38  km,  the  Sn 
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phase  is  very  prominent  at  the  shorter  ranges.  The  interference  of  reflected  and 
refracted  energy  gives  a  broad  locus  of  increased  amplitude  as  a  function  of 
distance.  However  once  the  velocity  jump  at  38  km  is  largely  eliminated  as  in 
model  SG3,  sharper  amplitude  maxima  occur  dominated  by  the  return  from  the 
crust-mantle  gradient  zone. 

The  amplitude  behaviour  of  these  theoretical  seismogram  calculations  as  a 
function  of  range  is  summarised  in  Figure  9.  The  maximum  amplitudes  are 
displayed  together  with  a  dashed  line  corres{X)nding  to  to  e"^  decay,  including 
the  effects  of  attenuation  in  the  model  at  the  2  Hz  centre  frequency  (7=  0.0028).  As 
would  be  expected  from  the  character  of  the  record  sections  in  Figure  8,  the 
calculated  amplitude  patterns  are  somewhat  erratic  and  cannot  be  fit  by  any  simple 
power  law  model.  The  model  SJO  with  a  sharp  Moho  transition  predicts  a  mean 
amplitude  decay  which  is  closest  to  r‘5/6  q-T  ^  ^ut  even  in  this  case  the  application 
of  a  simple  geometrical  spreading  correction  would  inevitaoly  impose  a  systematic 
pattern  of  amplitude  anomalies  on  the  corrected  data  being  used  to  estimate  Lg 
wave  attenuation  coefficients.  For  those  models  with  even  a  partial  gradient  zone 
the  rate  of  Lg  amplitude  decay  with  distance  is  enhanced  because  of  leakage  of 
energy  from  the  crustal  waveguide,  in  the  phase  velocity  range  contributing  to  Lg. 

An  additional  prediction  for  a  crust-mantle  transition  composed  almost  entirely 
of  a  gradient  zone  is  that  the  amplitude  maxima  for  Lg  are  sharply  peaked  as  a 
function  of  distance.  Model  SG3  predicts  very  str-^ng  Lg  amplitude  near  220  and 
440  km  and  it  is  interesting  to  note  that  in  the  normalized  seismograms  of  Figure  3 
the  Lg  arrivals  are  most  clear  at  just  these  distances.  The  crustal  models  of 
Finlayson  (1982)  designed  to  match  the  character  of  refraction  data  between 
Tennant  Creek  and  Mt.  Isa,  Queensland  show  lower  cnistal  velocities  exceeding 
7  km/s  and  gradients  rather  than  a  sharp  crust-mantle  transition.  Model  TCMl-3, 
in  particular,  has  an  average  gradient  (0.058  s"’)  between  35  and  60  km  which  is 
close  to  that  for  SG3  (0.064  s"’).  The  character  of  the  recorded  seismograms  along 
the  recorded  profile  and  the  observed  Lg  decay  rates  are  therefore  consistent  with 
the  pattern  expected  for  a  gradient  zone  between  crust  and  mantle. 

However,  our  comparison  of  different  crustal  structures  has  shown  the 
limitations  of  simple  corrections  for  geometric  spreading  in  estimating  attenuation 
coefficients  for  Lg.  Even  for  a  single  path  a  somewhat  erratic  trend  of  amplitude 
with  distance  can  be  expected,  albeit  with  a  superimposed  decay.  Such  decay 
behaviour  is  a  strong  function  of  the  character  of  the  crust-mantle  transition. 

Although  the  commonly  adopted  geometric  spreading  power  law  for  r  is 
seen  to  give  a  reasonable  representation  of  the  behaviour  for  a  sharp  Moho,  it 
cannot  be  universally  applied.  In  many  areas,  the  assumption  of  an  abrupt  crust- 
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mantle  transition  is  well  founded  but  beneath  Precambrian  terrains  (which  are 
well  represented  in  the  Soviet  Union)  a  gradational  crust-mantle  transition  seems 
to  be  quite  common. 

The  principal  application  of  Lg  attenuation  coefficients  lies  in  applying  path 
corrections  to  estimate  the  magnitude  of  a  source  (and  in  the  case  of  an  explosion 
to  thereby  estimate  the  yield).  The  use  of  inappropriate  amplitude  corrections 
could  lead  to  misleading  estimates  of  source  parameters. 

ACKNOWLEDGEMENTS 

We  thank  J.  Hulse  for  installing  the  portable  seismographs,  P.  Cummins  for 
programs  to  read  field  tapes,  S.  Ingate  and  K.  Muirhead  for  supplying  data  from  the 
Alice  Springs  array,  and  J.  Lahr  for  a  copy  of  HYPOELLIPSE.  This  work  was 
supported  in  part  by  the  Advanced  Research  Projects  Agency  of  the  U.S. 
Department  of  Defence  under  Grand  AFOSR-89-0330. 


REFERENCES 

Bolt,  B.A.  (1957).  Velocity  of  the  seismic  waves  Lg  and  Rg  across  Australia,  Nature, 
180, 495. 

Bowman,  J.R.  (1988).  Constraints  on  locations  of  large  intraplate  earthquakes  near 
Tennant  Creek,  Northern  Territory,  Australia  from  observations  at  the 
Warramunga  seismic  array,  Ceophys.  Res.  Lett.,  15,  1475-1478. 

Bowman,  J.R.,  G.  Gibson,  and  T.  Jones  (1990).  Aftershocks  of  the  22  January  1988 
Tennant  Creek,  Australia  intraplate  earthquakes:  Evidence  for  a  complex 
thrust  fault  geometry,  Ceophys.  /.  /Mf.,100,  87-97. 

Choy,  G.  and  J.R.  Bowman  (1990).  Rupture  process  of  a  multiple  main  shock 
sequence:  analysis  of  teleseismic,  local  and  field  observations  of  the  Tennant 
Creek,  Australia  earthquakes  of  22  January  1988,  in  press,  /.  Ceophys.  Res. 

Dwyer,  J.J.,  R.B.  Herrmann  and  O.W.  Nuttli  (1983).  Spatial  attenuation  of  Lg  waves 
in  the  central  United  States,  Bull.  Seism.  Soc.  Am.,  73,  781-796. 

Finlayson,  D.M.  (1981).  Reconnaissance  of  upper  crustal  seismic  velocities  in  the 
Tennant  Creek  block.  Northern  Territory,  BMR  J.  Aust.  Geol.  Ceophys.,  6, 
245-252. 

Finlayson,  D.M.  (1982).  Seismic  crustal  structure  of  the  Proterozoic  north 
Australian  craton  between  Tennant  Creek  and  Mount  Isa,  /.  Ceophys.  Res.,  87, 
10,569-10,578. 

Gaull,  B.A.,  P.J.  Gregson  and  K.  Malafant  (1989).  New  magnitude  scales  in  Westren 
Australia,  Aust.  Bur.  Min.  Resourc.  Record  1989/6. 

71 


Goleby,  B.R.,  R.D.  Shaw,  C.  Wright,  B.L.N.  Kennett  and  K.  Lambeck  (1989). 
Geophysical  evidence  for  "thick-skinned”  crustal  deformation  in  central 
Australia,  Nature,  337,  325-330. 

Jones,  T.D.,  G.  Gibson,  K.  McCue,  D.  Denham,  P.J.  Gregson,  and  J.R.  Bowman 
(1990).  Three  large  earthquakes  rupture  the  Australian  Precambrian  shield 
near  Tennant  Creek,  Northern  Territory  on  22  January  1988,  Geophys.  Res. 
Lett.,  in  press. 

Kennett,  B.L.N.  (1983).  Seismic  Wave  Propagation  in  Stratified  Media.  Cambridge 
University  Press,  Cambridge,  pp.  342. 

Kennett,  B.L.N.  (1988).  Systematic  approximations  to  the  seismic  wavefield,  in 
Seismological  Algorithms,  edited  by  D.J.  Doornbos,  Academic  Press,  London, 
237-259. 

Lahr,  J.C.  (1980).  HYPOELLIPSE/MULTICS;  A  computer  program  for  determining 
local  earthquake  hypocentral  parameters,  magnitude,  and  first  motion 
pattern,  U.S.  Geological  Survey  Open-File  Report  80-59,  31  p. 

Nuttli,  O.W.  (1973).  Seismic  wave  attenuation  and  magnitude  relations  for  eastern 
North  America,  /.  Geophys.  Res.,  78,  876-885. 

Plumb,  K.A.  (1979).  The  tectonic  evolution  of  Australia,  Earth  Sci.  Rev.,  14,  205-249. 

Pomeroy,  P.W.,  W.J.  Best  and  T.V  McEvilly  (1982).  Test  ban  treaty  verification 
with  regional  data-a  review.  Bull.  Seism.  Soc.  Am.,  72,  S89-S129. 

Street,  R.L.  (1978).  A  note  on  the  horizontal  to  vertical  Lg  wave-amplitude  ratio  in 
eastern  United  States,  Ea  thquake  Notes,  49,  15-20. 


72 


TABLE  1.  Earthquakes  used  in  this  study 


1  1 


No 

Origin  time 

latitude 

(OS) 

longitude 

(OE) 

depth 

(km)* 

Ml 

1 

3/27/88 

8  51  45.1 

19.853 

133.933 

4.0 

2 

3/27/88 

12  51  1.2 

19.785 

133.894 

4.0 

3.2 

3 

3/28/88 

0  27  12.6 

19.908 

134.030 

6.7 

5.0  (mb  4.9) 

4 

3/28/88 

0  58  14.0 

19.913 

134.028 

4.0 

3.5 

5 

3/28/88 

2  6  35.8 

19.795 

133.878 

4.0 

6 

3/28/88 

19  1  10.7 

19.858 

133.934 

4.0 

3.0 

7 

3/29/88 

13  22  40  3 

19.824 

133.923 

4.0 

8 

3/30/88 

15  54  39.5 

19.880 

134.004 

4.0 

3.1 

9 

3/31/88 

7  33  14.5 

19.827 

133.980 

4.0 

10 

3/31/88 

10  46  52.6 

19.899 

134.046 

4.0 

11 

4/  1/88 

5  9  58.1 

19.864 

134.101 

4.0 

*  constrained  to  4  km  unless  additional  local  data  were  available 


TABLE  2.  Crustal  velocity  model* 


Depth  (km) 

Vp  (km/s) 

0.00 

4.50 

0.17 

5.00 

0.34 

5.50 

2.14 

6.06 

6.00 

6.20 

13.00 

6.27 

27.00 

6.85 

40.00 

7.30 

47.40 

7.40 

50.50 

8.16 

58.00 

8.18 

61.00 

8.29 

*  Modified  from  Finlayson  (1981,  1982) 
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FIGURE  CAPTIONS 

1.  Map  of  the  source  (diamonds)  and  receiver  (triangles)  locations  used  for 
studying  Lg  propagation.  Portable  stations  (N1  and  H1-H8)  recorded  the  vertical 
component  only  at  25.6  samples/s  and  the  Alice  Springs  array  recorded  three- 
component  data  at  20  samples/s.  Stations  of  the  Warramunga  (WRA)  seismic 
array  are  shown  as  plus  signs.  Station  B2  at  the  intersection  of  the  two  arms  of 
the  WRA  array  recorded  three-component  data  at  16  samples/s. 

2.  Locations  of  the  aftershocks  used  for  the  Lg  stud/  A^ith  symbols  corresponding  to 
event  numbers  in  Table  1.  The  thrust  fault  scarps  (Bowman  et  al.,  1990)  are 
shown  as  heavy  lines  with  hatch  marks  on  the  upper  plate. 

3.  Seismograms  for  event  5  at  stations  on  the  line  between  Tennant  Creek  and 
Alice  Springs.  Each  trace  is  normalized  to  its  maximum  amplitude  which  is 
shown  above  in  digital  counts  for  stations  for  which  it  is  reliable.  Locations  of 
stations  are  shown  in  Figure  1,  and  station  SZ  is  the  short-period,  vertical- 
component  of  station  ASPA.  The  Lg  arrival  is  the  largest  phase  at  regional 
distances.  . 

4.  Three-component,  short-period  seismograms  for  event  5  from  the  Alice  Springs 
array.  Each  trace  is  normalized  to  its  maximum  amplitude,  which  is  shown  on 
the  right  in  digital  counts. 

5.  RMS  amplitude  for  2  Hz  center  frequency  and  corrected  for  r  geometric 
spreading  are  shown  as  a  funtion  of  source-receiver  distance  and  labeled  by 
event  number  as  in  Figure  2.  Lines  show  the  consistency  of  the  least  squares  fit 
to  the  data. 

6.  Attenuatior.  coefficient  y  as  a  function  of  the  power  n  assumed  for  geometrical 
spreading  for  1  Hz  (circles),  2  Hz  (triangles)  and  4  Hz  (diamonds).  RMS 
amplitudes  in  the  3.7-3. 1  km/s  group  velocity  window  were  used. 

7.  Velocity  models  used  for  synthetic  seismogram  calculations.  Heavy  solid  lines 
show  reference  model  SJO.  a)  Moho  discontinuity  replaced  by  a  gradient  in  SCil, 
SG2  and  SG3.  b)  Moho  discontinuity  reduced  in  models  SJ1,  SJ2  and  SG3. 

8.  Vertical-component  synthetic  seismograms  calculated  with  teh  wavenumber 
integral  method  (Kennett,  1983,  1988)  for  frequencies  1-  4  Hz,  with  explosive 
source  at  15  km  and  a  constant  Q  of  500  for  S  waves.  Reduction  velocity  is 
4  km/s.  a)  Gradient  models,  b)  Moho  jump  models. 

9.  Maximum  amplitude  of  synthetic  seismograms  shown  in  Figure  8  for  the 
velocity  models  in  Figure  7  as  a  function  of  source-receiver  distance.  Dashed 
lines  show  predicted  behaviour  for  an  Airy  phase  and  heavy  solid  lines  show 
amplitudes  for  model  SJO.  a)  Amplitudes  for  gradient  models  SGI,  SG2  and 
SG3.  b)  Amplitudes  for  jump  models  SJl  and  S]2. 
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